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ABSTRACT 

We study the effect of starliglit from tlie first stars on tlie ability of other 
minihaloes in their neighbourhood to form additional stars. The first stars 
in the ACDM universe are believed to have formed in minihaloes of total 
mass ^ 1O^~^M0 at redshifts z>20, when molecular hydrogen (H2) formed 
and cooled the dense gas at their centres, leading to gravitational collapse. 
Simulations suggest that the Population III (Pop III) stars thus formed were 
massive {'^ 100 A/q) and luminous enough in ionizing radiation to cause an 
ionization front (I-front) to sweep outward, through their host minihalo and 
beyond, into the intergalactic medium. Our previous work suggested that this 
I-front was trapped when it encountered other, nearby minihaloes, and that it 
failed to penetrate the dense gas at their centres within the lifetime of the Pop 
III stars ( < 3Myrs). The question of what the dynamical consequences were 
for these target minihaloes, of their exposure to the ionizing and dissociating 
starlight from the Pop III star requires further study, however. Towards this 
end, we have performed a series of detailed, ID, radiation-hydrodynamical 
simulations to answer the question of whether star formation in these sur- 
rounding minihaloes was triggered or suppressed by radiation from the first 
stars. We have varied the distance to the source (and, hence, the fiux) and 
the mass and evolutionary stage of the target haloes to quantify this effect. 
We find: (1) trapping of the I-front and its transformation from R-type to 
D-type, preceded by a shock front; (2) photoevaporation of the ionized gas 
(i.e. all gas originally located outside the trapping radius); (3) formation of an 
H2 precursor shell which leads the I-front, stimulated by partial photoioniza- 
tion; and (4) the shock- induced formation of H2 in the minihalo neutral core 
when the shock speeds up and partially ionizes the gas. The fate of the neu- 
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tral core is mostly determined by the response of the core to this shock front, 
which leads to molecular cooling and collapse that, when compared to the 
same halo without external radiation, is either: (a) expedited, (b) delayed, (c) 
unaltered, or (d) reversed or prevented, depending upon the flux (i.e. distance 
to the source) and the halo mass and evolutionary stage. When collapse is 
expedited, star formation in neighbouring minihaloes or in merging subhaloes 
within the host minihalo sometimes occurs within the lifetime of the first star. 
Roughly speaking, most haloes that were destined to cool, collapse, and form 
stars in the absence of external radiation are found to do so even when ex- 
posed to the first Pop 111 star in their neighbourhood, while those that would 
not have done so are still not able to. A widely held view that the first Pop 
III stars must exert cither positive or negative feedback on the formation of 
the stars in neighbouring minihaloes should, therefore, be revisited. 

Key words: cosmology: large-scale structure of universe - cosmology: theory 
— early universe - stars: formation - galaxies: formation 



1 INTRODUCTION 



Cosmological minihaloes at high redshift - i.e. dark-matter dominated haloes with virial tem- 
peratures Tvir < 10^ K, with masses above the Jeans mass in the intergalactic medium (IGM) 
before reionization (10^ < M/Mq < 10^) - are believed to have been the sites of the first star 
formation in the universe. To form a star, the gas inside these haloes must first have cooled 
radiatively and compressed, so that the baryonic component could become self-gravitating 
and gravitational collapse could ensue. For the neutral gas of H and He at T < 10'' K in- 
side minihaloes, this requires that a sufficient trace abundance of H2 molecules formed to 
cool the gas by atomic coUisional excitation of the rotational-vibrational lines of H2 . The 
formation of this trace abundance of H2 proceeds via the creation of intermediaries, or 
which act as catalysts, which in turn requires th e presence of a trace iq r iized fraction, in 
the following two-step gas-phase reactions (see, e.g., 



1967 



Lepp fc ShuU 



"SGB94"; 
H + e" 



1984 



Galh fc Palla 



H-+7, 



19981): 



Peebles fc Dicke 



1968; 



Shapiro fc Kang||l987l : IShapiro. Giroux. &: Babul 



Saslaw fc Zipoy 



I994J . henceforth, 
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(1) 



and 



H + H+^H+ + 7, 

H+ + H ^ H2 + H+. (2) 

Unless there is a strong destruction mechanism for H~ (e.g. cosmic microwave background 
at 2 > 100), the former (equation [1]) is generally the dominant process for H2 formation. 

Gas-dynamical simulations of the Cold Dark Matter (CDM) universe suggest that the 
first stars formed in this way when the dense gas at the centres of r ninihaloes of mass M ~ 
10^~^ Mfv, cooled and collapsed gravitationally at redshifts z ^ 20 fe.g. Abel. Brvan. fc Norman 



2000 



2001 



2002 



2003 



3romm. Coppi. fc 



Yoshida et al. 



iarson 



1999 



2002 



Yoshida et al. 



2003 



Machacek. Bryan, fc Abel 



20061 ). This work and others further suggest that these stars were 
massive (M* > 100 Mq), hot (Tcs — 10^ K), and short-lived (t* < 3 Myrs), thus copious emit- 
ters of ionizing and dissociating radiation. 

These stars constitute the Population III (Pop III) stars, or zero metallicity stars, which 
are believed to have exerted a strong, radiative feedback on their environment. The details 
of this feedback and even the overall sign (i.e. negative or positive) are poorly understood. 
Once the ionizing radiation escaped from its halo of origin, it created H II regions in the 
IGM, beginning the process of cosmic reionization. The photoheating which accompanies this 
photoionization raises the gas pressure in the IGM, thereby preventing baryons from collaps- 
ing gravitationally out of the IGM into new minihaloes when they form inside the H II re 



gions , an effect known as "Jeans-mass filtering" (SGB94: iGnedin fc Hui 



1998 



Oh fc Haiman 



20031 ). Inside the H II regions, whenever the I- 



minihaloes are subject to ph otoevaporation (jShapiro. Iliev. fc Raga 



ronts encounter pre-existing 



minihaloes, those 



2004, henceforth, SIR; 



Iliev. Shapiro, fc Ragall2005l . henceforth, ISR). A strong background of UV photons in the 
Lyman- Werner (LW) bands of II2 also builds up which can dissociate molecular hydrogen 
inside minihaloes even in the n eutral regions of the IGM, thereby disabling further collapse 
and, thence, star formation (e.g. 



Omukai fc Nishi 



1999 



Haiman. Abel, fc Rees 



200G 



Omukai 



20011 ). This conclusion changes, however, if some additional sources of partial ionization ex- 
isted to stimulate H2 formation without heating the gas to the usually high temperature of 
fully photoio nized gas (~ 10"^ K) at w'hich co Uisional dissociation occurs, such as X-rays from 



lar sources create a partially- ionized 



miniquasars (iHaiman. Rees. fc Loebl Il996af ) or if ste^ 
boundary layer outside of intergalactic H II regions (IRicotti. Gnedin. fc Shullll200ll ). Such 
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positive feedback effects, however, may have been only temporar y, because photoheating 



woul d soon become effective as background flux builds up over time ( iMacIntyre. Santoro. fc Thomas 
20M). 



The study of feedback effects has been limited mainly by technical difficulties. iHaiman et al. 



(I2OOOI ) studied the feedback of LW, ultraviolet (UV), and X-ray 



haloes without allowing hydrodynamic evolution. 



Ricotti et al. 



Dackgrounds on mini- 



feedback effect of stellar sources only on a static, uniform IGM. 



(12001) studied the radiative 



fl2002a 



Ricotti. Gnedin. &: ShuU 



b|) studied stellar feedback more self-consistently by performing cosmological hydro- 



dynamic simulations with radiative transfer, but the resolution of t h ese si mulations is not 



adequate for resolving minihaloes. iMachacek. Bryan. &: Abell (I2OOII . l2003l ) also performed 



cosmological hydrodynamic simulations, with higher resolution, but radiative feedback was 
treated assuming the optically thin limit, which overestimates the ionization efficiency, es- 
pecially in the high density regions which would initially be easily protected from ionizing 
radiation due to their high optical depth. The first self-consistent, radiation- hydrodynamical 
simulations of the feedback effect of external starlight on cosmological minihaloes were those 
of SIR and ISR, who studied the encounter between the intergalactic I-fronts that reionized 
the universe and individual minihaloes along their path. These simulations used Eulerian, 
grid-based hydrodynamics with radiative transfer and adaptive mesh refinement (AMR) to 
"zoom-in" with very high resolution, to demonstrate that the I-fronts from external ionizing 
sources are trapped when they encounter minihaloes, slowing down and transforming from 
weak, R-type to D-type, preceded by a shock. The gas on the ionized side of these I-fronts 
was found to be evaporated in a supersonic wind, and, if the radiative source continued 
to shine for a long enough time, the I-front eventually penetrated the minihaloes entirely 
and expelled all of the gas. These simulations elucidated the impact of the I-front and the 
physical effects of ionizing radiation on minihalo gas, quantifying the timescales and photon 
consumption required to complete the photoevaporation. They did not, however, address 
the aftermath of "interrupted" evaporation, when the source turns off before evaporation is 
finished. 



Recent studies by 



O'Shea et al. 



(l2005l ). lAlvarez. Bromm. fc Shapiro! (l2006al ). and lMesinger. Bryan. &: H 



(120061 ) addressed this question for minihaloes exposed to the radiation from the first Pop III 



star in their neighbourhood, instead of the effect of either a stea dily- driven I-front d uring 

( 20051 ) and 



global reionization or a uniform global background. The results of 



O'Shea et al. 
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Mesinger et al.l (120061 ) are seriously misleading, however, since they did not account properly 
for the optical depth to hydrogen ionizing photons. 



O'Shea et al 



(120051 ) assumed that the UV radiation from the first Pop III star that 
formed inside a minihalo in some region would fully ionize the gas in the neighbouring 
minihaloes. Using 3D hydrodynamics simulations, they found that, when the star turned 
off, H2 molecules formed in the dense gas that remained at the centre of the neighbouring 
minihalo, fast enough to cool the gas radiatively and cause gravitational collapse leading 
to more star for r nation . The H2 formation mechanism was the same as that described by 
Shapiro fc KangI (119871 ). in which ionized gas of primordial composition at a temperature 
T > 10^ K cools radiatively and recombines out of ionization equilibrium, enabling an en- 
hanced residual ionized fraction to drive reaction (1) (and [2], as well) as the temperature 



falls below the leve^ 
result 



O'Shea et al 



at which collisional dissociation suppresses molecule formation. As a 
1 

(120051 ) concluded that the radiative feedback of the first Pop III stars 
was positive, triggering a second generation of star formation in the minihaloes surrounding 
the one that hosted the first star. 



Mesinger. Bryan, fc HaimanI (l2006l ) also used 3D hydrodynamics simulations to consider 
the fate of the gas in the relic H II regions created by the first Pop III stars, they concluded 
that the radiative feedback of the first stars could be either negative or positive and estimated 
a critical UV intensity which would mark the transition from negative to positive feedback. 



Mesinger et a 



been do ne by 
effect in 



(|2006r). howeve r , stud i ed thi s effect only in the optically thin limit, as had also 



Machace 



O'Shea et al. 



c et al 



3) 



(12001 



(120051 ) and 



2003!). The main m echanisms of the positive feedback 



Mesinger et al.l (120061 ) are, therefore, identical. 



Alvarez et al. 



(l2006al ). on the other hand, performed a high-resolution ray-tracing calcu- 
lation to track the position of the I-front created by the first Pop III star as it swept outward 
in the density field of a 3D cosmological SPH simulation of primordial star formation in the 
ACDM universe over the lifetime of the star. When this I-front encountered the minihaloes 
in the neighbourhood of the one which hosted the first Pop III star, it was trapped by the 
minihalo gas before it could reach the high-density region (core), due to the minihalo's high 
column density of neutral h ydrogen. This i s consi stent with the results of SIR and ISR men- 



tioned above. According to 



Alvarez et al. 



(l2006al ). in fact, the lifetime of the Pop III star is 
less than the evaporation times determined by SIR and ISR for the relevant minihalo masses 
and flux levels in this case, so the neutral gas in the core is never ionized by the I-front. It 
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seems that the initial assumption of full ionization of nearby h aloes by 



and the optically thin limit assumed by iMesinger et al.l (120061 ) are invalid. 



O'Shea et al. 



(120051) 



The final fate of this protected neutral cor e, however, is still unclear, because the I-front 



tracking calculations by 



Alvarez et al. 



(j2006af ) did not include the hydrodynamical response 
of the minihalo gas to its ionization, a full treatment of radiative transfer or the primordial 
chemistry involving H2. One might naively expect that the nett effect would be negative, 
because heating from photoionization would ultimately expel most of gas from minihaloes, 
although the results of SIR and ISR, again, show that this minihalo evaporation would not be 
complete within the lifetime of the Pop III star. On the other hand, partial ionization beyond 
the I-front by hard photons from a Pop III star might be able to promote II2 formation, 
once the dissociating UV radiation from the star is turned off, which would then lead to a 
cooling and collapsing core. This issue can be addressed only by a fully coupled calculation 
of radiative transfer, chemistry, and hydrodynamics, which will be the focus of this paper. 

We shall attempt to answer the following questions: Does the light from the first Pop III 
star in some neighbourhood promote or prevent the formation of more Pop III stars in the 
surrounding minihaloes? More specifically, do the neutral cores of these nearby minihaloes, 
which are shielded from the ionizing radiation from the external Pop III star, subsequently 
cool and collapse gravitationally, as they must in order to form stars, or are they prevented 
from doing so? Towards this end, we simulate the evolution of these target haloes under 
the influence of an external Pop III star using the 1-D spherical, Lagrangian, radiation- 
hydrodynamics code we have developed. We adopt a 120 Mq Pop III star as a source, and 
place different mass haloes at different distances to explore a wide range of the parameter 
space for this problem. Masses of target haloes are chosen to span the range from those too 
low for haloes to cool and collapse by H2 cooling without external radiation to those massive 
enough to do so on their own. 

Our calculation is the first self-consistent gas-dynamical calculation of the feedback 
effects of a single Pop III star on nearby haloe s. A similar approach by 1-D radiation- 



Kitayama et al 



(I2OOII ) . Their work, how- 



hydrodynamics calculation has been performed by 
ever, focuses on the effect of a steady global background from quasars and from stars with sur- 
face temperatures T^, ~ 10^ K, rather than a single, short-lived Pop III star with ~ 10^ K. 



In addition, whi 
was reported by 



e we were prepari n g this manuscript, a study which is similar to our work 



Susa fc Umemural (120061 ) . where a 3D radiation- hydrodynamics calculation 
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with SPH particles was performecll]. A major difference of their work from ours is that they 
focus on the subclumps of the halo which hosts the first Pop III star, while we focus on ex- 
ternal minihaloes in the neighbourhood of such a host halo. We also apply a more accurate 
treatment of H2 self-shielding, as well as a more complete chemistry network of neutral and 
ionic species of H, He, and H2. A more fundamental difference from these previous studies is 
our finding of a novel H2 formation mechanism: collisional ionization of pre- 1- front gas by a 
shock detached from a D-type I-front. This mechanism occurs at the centre of target haloes, 
which would otherwise remain very neutral. This mechanism creates new electrons abundant 
enough to promote further H2 formation, which can even expedite the core collapse. 

In Section [2] we describe the details of the 1-D spherical radiation-hydrodynamics code 
we have developed. Some details left out in Section [2] will be described in Appendices. In 
Section 13.11 we describe the initial setup of our problem. We briefly describe a test case in 



Section m where we let a m inihalo evolve 



rom an initially ionized state, to show that our 



code reproduces the result of lO 'Shea et al.l (120051 ) in that case. In Section [5] and Section [6l we 



present the main results of our full radiation-hydrodynamics calculation. We summarize our 
results in Section[7l Throughout this paper, we use the ACDM cosmological par ameters, (Vt\, 



h) = (0.73, 0.27, 0.043, 0.7), consistent with the ITM^P first-year data (jSpergel et al. 

2003k 



2 NUMERICAL METHOD: 1-D SPHERICAL, 

RADIATION-HYDRODYNAMICS WITH PRIMORDIAL CHEMISTRY 
NETWORK 

In this section, we describe in detail the 1-D spherical, Lagrangian, radiation-hydrodynamics 
code we have developed for both dark and baryonic matter. We describe how hydrodynamics, 
dark matter dynamics, radiative transfer, radiative heating and cooling, and finally the 
nonequilibrium chemistry are handled. The finite differencing scheme, reaction rates, and 
certain other details not treated in this section will be described in Appendices. We include 



^ A new preprint bv lAbel. Wise, fc BrvarJ 1 20061 ) has also appeared which addresses this issue. We will discuss this further in 
Section 16.61 

^ As we do not perform a statistical study, our result is independent of the cosmic density power spectrum. The three-year 
WMA P data does not show a big discrepancy in the set of cosmological parameters of the interest in this paper i Spergel et al. 
20061 ). The change in trg an d the index of the pri mordial power spectrum n would translate to ^ 1.4 redshift delay of structure 



formation and reionization ( Alvarez et al 



2006d) 
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the neutral and ionic species of H, He and H2, namely H, H"*", He, He"*", He"'""'", H~, H2, Ha" 
and e~, in order to treat the primordial chemistry fully. As deuterium and lithium exist in 
a negligible amount, we neglect D and Li specieqf . 



2.1 Hydrodynamic Conservation Equations 

The baryonic gas obeys inviscid fluid conservation equations, 

+ + ) + = -P-T^ (4) 

De p d ^ T-A 
Dt pr^or p 

where e = (3p) / (2p) is the internal energy per unit baryon mass, F is the external heating 
rate, and A is the radiative cooling rate. Note that all the variables in equations ([3]) - ([5]) 
denote baryonic properties, except for m, the mass enclosed by a radius r, which is composed 
of both dark and baryonic matter. 

We do not change the adiabatic index 7 throughout the simulation. As long as monatomic 
species, H and He, dominate the abundance, 7 = 5/3 is the right value to use. This ratio of 
specific heats, 7, can change significantly, however, if a large fraction of H is converted into 
molecules. For example, the three-body H2 formation process, 

H + H + H ^ H2 + H, (6) 

will occur vigorously when nn ^ 10^ cm s~^ and T < 10^ K, which will invalidate the use of a 
constant 7. To circumvent such a problem, when such high density occurs, we simply stop 
the simula tion. This proces s is, nevertheless, important in forming the protostellar molecular 
cloud (e.g. 



AbeletaL 



2OO2I ) . This issue will be further discussed in Section |6l when we define 



the criterion for the collapse of cooling regions. 



D and Li components have usually been neglected du e to their relatively low abundance, hence the n egli gible contribution to 



cooling (e.g. lLepp h Shulllll984l : IShaDiro fc Kan^[l987l) . Recent stud es by 



Nagakura &: Omukai 



|2005j) and 



Johnson fc Bromm 



I 2006), however, show that enough HD is generated in strongly-shocked, ionized primordial gas which then can cool below 
the temperature of ~ 100 K already achieved b y H2 cooling alone, down t o the temperature of the CMB. As the HD cooling 
process is negligible if gas remains neutral (e.g. Ijohnson &: BromrrJ 1 20061 )). however, we may neglect the HD cooling process 
in our calculation as long as we are interested in the centre of target haloes which remains mostly neutral at any time. We will 
discuss this issue further in Section [7] 
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The shock is treated using the usual artificial viscosity technique (e.g. 



Von Neumann fc Richtmyer 



1950). The pressure p in equations (jl]) and ([5]) contains the artificial viscosity term. The de- 
tails of this implementation are described in Appendix A. 



2.2 Dark Matter Dynamics 

Gravity is contributed both by the dark matter and the baryonic components. Let us first 
focus on the dark matter component. In order to treat the dark matter gravity under spher- 
ical symmetry, almost all previous studies have used either a frozen da rk matter potential 



Thoul &: Weinber 



or a set of self-gravitating dark matter shells in radial motion only (e.g. 
19951 ). Both methods have their own limitations. The frozen potential approximation can- 
not address the effect of a possible evolution of the gravitational potential. The radial-only 
dark matter approximation suffers from the lack of any tangential motion, producing a viri- 
alized structure whose central density profile is much steeper (p oc r~^ with /3 > 2; see 
e.g. ) than that of haloes in c osmological, 3-D N-body simulations (/? ~ 1, as found in 
Navarro. Frenk. fc Whitel 119971 ). 

In order to treat the dynamics of dark matter more accurately than these previous 
treatments, we use th e the fiuid approximation we have developed and reported elsewhere 



(lAhn fc Shapiro 



see 



2005f). W e briefiy summarize its derivation here; for a detailed description. 



Ahn fc Shapiro! (120051 ). Collisionless CDM particles are described by the collisionless 



Boltzmann equation. When integrated, it y ields an infinite set of con servation equations, 
which is called the BBGKY hierarchy (e.g. iBinney fc Tremaind 119871 ). However, CDM N- 
body simulations show that virialized haloes are well approximated by spherical symmetry. 
These simulations also show that the velocity dispersions are highly isotropic: radial disper- 
sion is almost the same as the tangential dispersion. These two conditions make it possible 
to truncate the hierarchy of equations to a good approximation, which then yields only 
three sets of conservation equations. Amazingly enough, these equations are identical to the 
normal fiuid conservation equations for the adiabatic index 7 = 5/3 gas: 

0, (7) 

Gm 



9pd . d 2 



d d 2 

-rrXpdUd) + ^ij>d + Pdul) + -pdul 
ot or r 



-pd- 



Pd d 



(8) 
(9) 



Dt pd r'^dr 

where the subscript d represents dark matter, the effective pressure pd = Pd {ud ~ {ud))'^ is 
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the product of the dark matter density and the velocity dispersion at a given radius, and 
the effective internal energy per dark matter mass = 3pd/2pd. We use these effective fluid 
conservation equations (equation d [HI EI) to handle the motion of dark matter particles. 

Note that dark matter shells in this code represent a collection of dark matter particles in 
spherical bins, in order to describe "coarse-grained" properties such as density (p^) and the 
effective pressure (p^). As these coarse-grained variables follow the usual fluid conservation 
equations, the hyperbolicity of these equations leads to the formation of an effective "shock." 
The location of this shock will determine the effective "post-shock" region. This post-shock 
region corresponds to the dark matter shell- crossing region. Because of the presence of this 
effective shock, we also use the artificial viscosity technique. This collisional behaviour of 
our coarse-grained dark matter shells originates from o ur ch oice of physical variable. For 



Ahn fc Shapirol (120051 ) and lAlvarez et al.l (120031 ) for 



further details, the reader is referred to 
description and application of our fluid approximation. 
The mass enclosed by a dark matter shell of radius r, 

m(< r) = mDM(< r) + mbary(< r), (10) 

enters equations (jll) and (IHl). When computing m(< r), we properly take account of the 
mismatch of the location of dark matter shells and baryon shells. 



2.3 Radiative transfer 

A full, multi-frequency, radiative transfer calculation is performed in the code. Since H2 
cooling is of prime importance here, we first pay special attention to calculating the optical 
depth to UV dissociating photons in the LW bands and the corresponding H2 self-shielding 
function. We then describe how we calculate the optical depth associated with any other 
species depending upon the location of the radiation source. The finite difference scheme for 
the calculation of radiative rates is described in the Appendix A. 

2. 3. 1 Photodissociation of H2 and Self- Shielding 

Hydrogen molecules are photodissociated when a UV photon in the LW bands between 11 eV 
and 13.6 eV excites H2 to an excited electronic state from which dissociation sometimes 
occurs. When the column density of H2 becomes high enough (A^H2 ^ 10^'' cm~^), the optical 
depth to photons in these Lyman- Werner bands can be high, so H2 can "self-shield" from 
dissociating photons. Exact calculation of this self-shielding requires a full treatment of all 76 
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Lyman- Werner lines, even when only the lowest energy level transitions are included. Such 
a calculation is feasible under simplified conditions such as a radiati ve transfer problem 



through a static medium (e.g. 



Haiman et al 



2000 



Ricotti et al. 



2OO1I ). Unfortunately, for 



combined calculations of radiative transfer and hydrodynamics, such a full treatment is 
computationally very expensive. 

Under certain circumstances, however, one can use a pre-computed self-shielding function 
expressed in terms of the molecular column density A/^Hj and the temperature T of gas, which 
saves a great amount of computation ti me. In a cold, static medi um, for instance, one can 



use a self-shielding function provided by 



Draine fc Bertoldil fll996h : 



Ffih\ 



shield 



mm 



,10i4cm-2. 

The photodissociation rate is then given by 



ku^ = 1.38 X 10^ (Ju) 



hu=12.87cY -^shield, 



(11) 



(12) 



where (J^^ 



hu=12.87cY 



(ergs ^ cm ^ Hz ^ sr ^) is the mean intensity in the spectral region 



of the LW bands. This a pproximation has bee i i widely used in 



structure forma 


;ion ( 


Kitavama et al. 


2004) 



e.g. 



Kitayama et al. 



2001 



Glover fc Brand 



;he s t udy of high redshift 



2001 



Yoshida et al. 



2003; 



The problem with equation flTTl) is that when the gas temperature is high or gas has 
motion along the line of sight to the source, the thermal and velocity broadening of the 
LW bands caused by the Doppler effect can signifi cantly reduce the optical d epth. A better 



Draine fc Bertoldil (119961 ) . now in terms 



treatment for thermal broadening is also given by 
of the molecular column density A^H2 and the velocity-spread parameter b of the gas: 
0.965 0.035 



Ffih\ 



shield 



1 



x/b5] 



+ 



2 (l+x)0-5 
4/1 I ™\0.5i 



X exp[-8.5 X 10-^(1 + xy-% (13) 

where x = Nnj5 x 10^^ cm'^, 65 = 6/ 10^cms~\ and b = 1.29 x 10^ (Tk/A)^^^ cms~^ where 
A is the atomic weight Jspitzeilll978l ). For H2, b = 9.12 km s'^ (T/IO^K)^/^ 

In the problem treated in this paper, we frequently find T ^ 10^ — 5 x 10'^ K in the 
gas parcel (shell) which contributes most of the H2 column density. We also find that this 
gas parcel usually moves at ~ 2 — 5kms~^ (see Section 16. 2p . The combined effect of 
the thermal broadening and the Doppler shift on the shielding function, then, may be well 
approximated by a thermally broadened shielding function with T ^ 10'^ K. Throughout 
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Figure 1. Power-law self-shielding function for cold, static gas vs. self-shielding function for hot gas at T = 10* K. The problem 
of interest to us resides in the sensitive region, A'^h, ~ 10^* — 10^® cm~^, where the biggest discrepancy exists. 

this paper, therefore, we use equation (fT3!) with T = 10^ K to calculate the self-shielding. 
For the photo-dissociation rate, we use equation f[T^ . 

We show in Fig. [T] how much the static, cold shielding function (equation [TT]1 may over- 
estimate the self-shielding in our problem, by comparing this to the thermally-broadened 
shielding function (equation [T3l) at T = 10^ K. The biggest discrepancy between these two 
shielding functions exists for A'^h2 ~ 10^^ — lO^^cm"^. Interestingly enough, the H2 column 
density in our problem usually resides in this regime. It is crucial, therefore, to take into 
account the effects of thermal broadening and Doppler shift carefully, as we do in this paper. 



2.3.2 External Source 

Since our calculations are 1-D, spherically-symmetric, we have assumed the external radia- 
tion source contributes a radial flux F°^*(r) at frequency u and radius r, measured from the 
minihalo centre, given by 



-T„{>r) 



(14) 



where L^^* is the source luminosity, and rj^(> r) is the optical depth along the radial direction 
from radius r to the source located at a distance r = D. 

The radiative rate of species i at radius r is then given by 

(15) 
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where we have used the fact that AnJ,^ = -F^^*, as long as the external radiation can be 
approximated as a ID planar flux. In practice, one calculates this rate in a given grid-cell - 
i.e. spherical shell - with finite thickness. If such a grid-cell has a small optical depth, F^^^ is 
almost constant across the grid, so one could take the grid-centred value of F^;^^ to calculate 
ki{r). This naive scheme, however, does not yield an accurate result when a grid-cell is 
optically thick, where F^^^* may vary significantly over the cell width. This problem occurs 
frequently for solving radiative transfer through optically thick media, where individual cells 
have large optical depth. In order to resolve thi s prob lem, we use a "photon -conserving" 



scheme like that described by 



Razoumov fc Scott 



(jl999|) and lAbel et all (119991 ). The details 



of our implementation of this scheme are described in the Appendix A. 

2.4 Heating and Cooling 

2.4- i Photoheating 

Photoheating results from thermalization of the residual kinetic energy of electrons after 
they are photoionized. In general, the photoheating function is described by 



oo 



i i ^ 

= Vnj / du^ — -{hv ~ huixY,), (16) 

where /ii^j^th is the threshold energy over which the residual photon energy is converted into 
the kinetic energy of electrons, and the nett heating function F is the sum of individual heat- 
ing functions ({Fj}). In finite-differencing equation f|T6|) . we also use the photon-conserving 
scheme as we do for equation (fT5l) . This prevents cells with large optical depth from obtaining 
unphysically high heating rates. See Appendix A for details. 

2.4-2 Radiative cooling 

Cooling occurs through various processes. For atomic species, it comes from coUisional exci- 
tation, coUisional ionization, recombination, free-free emission, and CMB photons scattering 
off free electrons (Compton cooling/heating). For atomic H and He, cooling is dominated 
by coUisional excitation (for T<2 x lO^K) and free-free emission (for T>2 x lO^K). The 
atomic cooling rate decreases rapidly at T < lO^K, as there are no collisions energetic enough 
to cause excitation. It is difficult, therefore, to cool gas below T lO^K solely by atomic 
cooling of primordial gas. 

© 0000 RAS, MNRAS 000, 000-000 



14 K. Ahn & P. R. Shapiro 

Molecular hydrogen (H2), however, is able to cool gas below T ^ lO^K, down to T ^ 
lOOK, by collisional excitation of rotational-vibrational lines by H atoms. An important 
question to address is how much H2 is created, maintained, or destroyed under the influence 
of an ionizing and dissociating radiation field. E ven a small fraction, riHo/nn 10~^, is 



sometimes enough to cool gas below IC^K (e.g. see IShapiro &: Kang 



We use cooling rates in the parametrized forms given by 



1987h 



Anninos et al. 



the hydrogen molecular cooling. For H2 cooling, we use the fit given by 



19971). ex c ept fo r 



Galh fc Fallal (119981 ) ■ 



where t he low density coolin g rate has been updated significantly from the previously used 



rate by lLepp &: Shulll (1l984l ). which suffers from the uncertainties associated with the only 
collisio nal coefficients available at that time. At low densities, rin ^ lO^cm"' ^, the cooling 



rate of 



Lepp &: Shulll (119841 ) is bigger by an order of magnitude than that of 



Galh fc Falla 



mm at T ^ lOOOii'. 



2.5 Nonequilibrium chemistry 



The general rate equation for the abundance of species i is given by 
Qyi ■ 

= Ci{T, {rij}) - Di{T, {nj})ni, 



(17) 



where Ci is the collective source term for the creation of species i, and the second term is the 
collective "sink" term for the destruction of species i. The processes included and adopted 
are shown in Table 



by 



Shapiro fc Kang 



Bll in Appendix B. Most of the rate coefficients are those from the fits 



(Il987l ). with a few updates. 



We also adopt the rate solving scheme proposed by 



Abel et al. 



(I1997I ). It is well known 
that coupled rate equations in the form of equation (fTTl) are "stiff" differential equations, 



Abel et al. 



whose numerical solution suffers from instability if explicit ODE solvers are used. 
(I1997I ) show that their implicit, backward difference scheme provides enough stability. Accu- 
racy of the solution is achieved by updating each species in some specific order, rather than 
updating all species simultaneously from their values at the last time step. In addition, the 
abundance of the relatively fast reactions of H~ and H.2 are approximated by their equilib- 
rium values, which are expressed by simple algebraic equations. See the Appendix A for the 
corresponding finite-differencing scheme. 



We will frequently quote our results in terms of the fractional number density of species 



where nn is the number density of the total atomic hydrogen atoms. We use 
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X, however, to denote the fractional electron number density, i/e, which is a measure of the 
ionized fraction. 



2.6 Code tests 



We tested our code against the following problems which have analytic solutions: 

(A) the self-similar, spherical, cosmological infall and accretion shock resulting from a 
point-mass perturb ation in an Einstein-de Sitter universe of gas and coUisionless dark matter 
( lBertschingeiill985l ): 



(B) the self-similar blast wave w 



lich results from a strong, adiabatic point explosion in 



a uniform gas - the Sedov solution (ISedovlll959l ) 



(C) the propagation of an I-front from a steady point-source in a uniform, static medium 

(D) the g as-dynamical expansion of an H II region from a point source in a uniform gas 



( iLasker 



19661 ) 



(E) the gas-dynamical expansion-phase of the H II region from a point-source in a 



nonuniform gas whose density varies with d istance r from the source as r 



( iFranco. Tenorio-Tagle. fc Bodenheimer 



199o[ ). 



w 



3/2 



Our code passed all the tests described above with an acceptable accuracy. Test results 
are described in Appendix C. 



3 THE SIMULATIONS 



3.1 Initial Setup 

We now describe the initial setup for the problem of radiative feedback effects of Pop III 
stars on nearby haloes at z ^ 20. The first stars form inside rare, high density peaks at 
high redshift. We place target haloes of different mass M = [2.5 x 10"^, 5 x 10^, 10^, 2 x 
10^, 4 X 10^, 8 X 10^] Mq at different locations from the source, with proper distance D = 
{180, 360, 540, 1000} pc, which are all assumed to be affected directly by the radiation 



field from the source Pop III star of mass = 120 Mq 



We ex pose the target halo to this 



radiation field for the lifetime of the star, t*(120 Mq) ~ 2.5 Myrs (jSchaerer 



2OO2I). The source 



^ The additional case of D = 50 pc, Fq = 600, M = 5.5 X 10^ Mq, will be discussed separately in Section 16.61 with regard to 
the case in which the target minihalo is merging with the minihalo which hosts the star, separated by less than its virial radius 
from the star 
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Pop III star is assumed to be located in a halo of mass M ~ IO^Mq. Time is measured from 
the arrival of the stellar radiation at the location of the target minihalo. 



This setup is well justified by the cosmological simulations by lAlvarez et al.l (l2006al ). A 
cosmological gas and N-body simulation of structure formation in the ACDM universe on 
small scales by a GADGET/SPH code was used to identify the site at which the first Pop 
III star would form. This occurred at 2; = 20, at the location of the highest density SPH 
particle in the simulation box, located within a halo of m ass M lO^M f r,. This provided the 



initial density field for the I-front tracking calculations in 



Alvarez et al. 



(I2nn6al ). The I-front 



from this first star escaped from the host halo quickly with high escape fraction, traveling 
as a supersonic, weak R-type front. By the end of the lifetime of the star (~ [3 — 2] Myrs) 
for stellar masses in the range M^, ~ [80 — 200] M©, the star's H II region had reached a 
maximum radius of about 3 kpc. 

We approximate the spectral energy distribution (SED) of the source star by a blackbody 
spectrum. A Pop III star of mass M^, ^ 120 Mq, according to ISchaererl (120021 ). has the time- 
average effective temperature T^g ^ lO^K and luminosity L = dvL^ ^ lO^'^^^L©. The cor- 
responding ionizing photon luminosity with this blackbody spectrum is = dvL^/hv = 
1.5 X 10^'^s^^, where hun = 13.6 eV is the hydrogen ionization threshold energy. We assume 
that the source radiates with these time-averaged values throughout its lifetime, then stops. 
As the photons escape in a time scale short compared to the lifetime of the star and the 
escape fraction is high, we simply ignore the effect of the intervening gas (e.g. optical depth 
from the host halo and the IGM) and assume that the bare radiation field hits the edge of 
target haloes directly. 

As we fix the luminosity of the source, different distances correspond to different fluxes. 
We express the frequency-integrated ionizing photon flux, F in units of 10^° s^^ kpc^^, to 
give the dimensionless flux, Fq = A^ph,5o/-D^pc = A^ph.seZ-DMpc; where iVph,5o is the ionizing 
photon luminosity (in units of 10^° s^^) and -Dkpc (-Dmpc) is the distance in units of kpc 
(Mpc), respectively. The value Fq ~ 1 is typi cal for minihaloes e ncountered by intergalactic 
I-fronts during global reionization (e.g. see IShapiro et al.ll2004j ). Interestingly enough, Fq 
for our "small-scale" problem has a similar value. The Pop III star in our problem has 
A^ph 50 = (5*/10^° = 1.5. For distances 180 pc, 360 pc, 540 pc and 1000 pc, Fq corresponds 
to 46.3, 11.6, 5.14 and 1.5, respectively. 
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3.2 Initial Halo Structure 

For the initia l halo structure, we adopt the minimum-energy tru ncated isothermal sphere 
(TIS) model (jshapiro. Iliev. fc Raga IQQol : Iliev fc Shapiro 2001 ). which will be described 



further in Section I3.2.1[ The thermodynamic properties and chemical abundances of the gas 
in these target haloes, however, is somewhat ambiguous. The density and virial temperature 
of these haloes are higher than those of the IGM in general, which drives their chemical 
abundances to change from the IGM equilibrium state to a new equilibrium state. The 
most notable feature is the change of ?/h2 ^iid x. The IGM equilibrium value of the electron 
abundance, x ~ 10~^, is high enough to promote H2 formation inside minihaloes to yield a 
high molecule fraction, ^ 10~^ — 10~^. At the density of gas in the halo core, this newly 
created H2 is capable of cooling the minihalo gas to T ^ 100 K, and depending on the virial 
temperature, the minihalo may, therefore, undergo a runaway collapse. 

The time for this evolution of the target halo gas is short compared to the age of the 
universe when the first star forms in their neighbourhood. As a result, it is likely that the 
target haloes are exposed to the ionizing and dissociating radiation from that first star as they 
are in the midst of evolving, with fine-tuning required to catch all of them in a particular stage 
of this evolution. As the evolutionary "phase" of our target haloes is uncertain, we adopt two 
different phases as our representative initial conditions. In Phase I, chemical abundances have 
not yet evolved away from their IGM equilibrium values. This stage is characterized by low 
H2 fraction, y^^ ~ 2 x 10^^ and high electron fraction, x ~ 10~^. Phase II is the state which 
is reached, after allowing the Phase I minihalo to evolve chemically, thermodynamically and 
hydrodynamically for a few million years (a small fraction of a Hubble time, = 186Myrs 
at z = 20), until the electron fraction has decreased to x ~ 10^^. Phase II is characterized 
by high H2 fraction, y^,^ ~ 10"'^ — 10^^, and cooling-induced compression of the core relative 
to Phase I, by a factor between 1 and 20, higher for higher minihalo mass. 



3.2.1 Phase I: Unevolved Halo with IGM chemical abundance in hydrostatic equilibrium 

The first phase we choose is the initial state we assumed above, namely the nonsingular TIS 
structure with IGM chemical abundances. This phase is characterized by gas in hydrostatic 
equilibrium, with the truncation radius (outer boundary of the halo) 

n = 102.3^ — 
V0.27/ \0.7J 
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M 



1/3 



\2 ■ IO^MqJ 
the virial temperature 



l + z 
1 + 20 



T 



593.5 



1.22 

M 







0.27 

2/3 



1/3 



0.7 



pc, 



2/3 



(18) 



X 



2 ■ lO^Mm 



l + z 
1 + 20 



K, 



(19) 



where fi is the mean molecular weight (1.22 for neutral gas and 0.59 for ionized gas) and 
the central density 



Po 



4.144 X 10" 



-22 



^^0 

0.27 



0.7, 



l + z 



g cm 



1 + 20. 

which can also be expressed in terms of the hydrogen number density by 

x{n,/no)po 



(20) 



30 



X 



0.76/ V0.043 



0.7, 



l + z 
1 + 20 



cm 



(21) 



where X is the hydrogen mass fraction in the baryon component. This central density is 
about 1.8 X 10^ p(z), where p{z) is the mean ni a tter d ensit y at redshift z, while a t r = Vtr, 
p = 35p(2;). For more details, see IShapiro et al.l (1l999l ) and llliev &: Shapiro! (1200 ll ). 

We assign chemical abundances that reflect the IGM equilibrium state, which is char- 
acterized by high electron fraction - high enough to promote H2 formation under the right 
conditions - and low H2 fraction - low enough to contribute negligible molecular cooling. 



We adopt yn = I, yue = 
species (see, e.g. SGB94; 



0.0789, X 



Ricotti et al. 



10-4, = 2 X 10"^ and {yi} = for other 



2OO1I). 



3.2.2 Phase II: Evolved Halo with Recombimng and Cooling Core 

The second initial condition we choose is the evolved state (Phase II) reached by allowing 
the system to evolve from Phase I initial conditions before the arrival of radiation from the 
Pop III star. In particular, we follow this evolution until the central electron fraction has 
dropped to 10^^ by recombination from Phase I. We choose this condition because it is now 
characterized by high molecule and low electron fraction, contrary to Phase I. The fate of 
this halo will then mainly be determined by how easily this abundant II2 is protected against 
dissociating radiation after the star turns on. The answer will also depend upon how much 
change has occurred hydrodynamically, because in some cases the halo core may have cooled 
and collapsed significantly enough to be unaffected by the feedback from late irradiation. 
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Figure 2. Initial conditions for target lialocs. We clioosc two difFcrcnt pliascs of TIS iialo evolution as separate initial conditions. 
Phase I (unevolved; IGM abundance; solid) and Phase 11 (evolved from Phase 1 for a time Atj n until x = 10~® at centre; 
dotted) are plotted for each mass of target halo. Each panel is labelled with the value of Ati^u in Myrs. Note that hydrodynamic 
diflference between two phases is evident in haloes of mass M >2 x 10^ Mq . 
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The time to reach Phase II is different for different mass haloes because of different 
gas properties. Initially, as we start from the TIS density profile whose central density is 
independent of the halo mass, the recombination rate is higher for smaller mass haloes, 
because hydrogen recombines according to the following: 

^ oc nun,-T-'-\ (22) 

The situation becomes complicated, however, once evolution begins and density changes. The 
H2 cooling and collapse in the central region of the haloes is increasingly effective as halo 
mass increases, because of the increasingly large difference between the virial temperature 
and the II2 cooling temperature plateau, ~ 100 K. The corresponding rapid collapse and 
cooling in massive haloes can easily offset the initial temperature dependence by obtaining 
high density and low temperature, as is seen in equation (l22l) . Phase II for large mass haloes 
represents haloes that have already started their cooling and collapse. 

In Fig. [21 we show halo profiles in Phase I and Phase II for different halo masses. We also 
show how much time it takes for the haloes to evolve from Phase I to Phase II. The times 
for gas at the halo centre to recombine to x = 10~^ are in the range 7 < Ati ii(Myrs) < 24 
for halo masses 0.25 < M/{10^Mq) < 8, peaked at Ati,ii = 24Myrs for M = 5 x W Mq. 
In all cases, Atj^n <^ = 186Myrs, the age of the universe at z = 20. 



4 HALO EVOLUTION FROM FULLY-IONIZED INITIAL CONDITIONS: 
THE CONSEQUENCES OF IRRADIATION WITHOUT OPTICAL 
DEPTH 

Before describing the results of our full radiative transfer, hydrodynamics calculation, we 
describe an experiment designed to show the effect of neglecting the optical depth of the 
minihalo to ionizing radiation from the external star during the star's lifetime on the mini- 
halo's evolution after the star shuts off. For this purpose, we assume the target minihalo is 
initially fully-ionized and heated to the temperature of a photoionized gas as it would be if 
it were instantaneously fla sh-ionized by starlig ht in the optically-thin limit. Such a setup is 



equivalent to that used by 



O'Shea et al 



(I2OO5I ). where they find that second-generation star 



formation is triggered when the ionization of the minihalo caused by the nearby Pop III star 
leads to cooling by H2. The high initial electron fraction is present because of the assumption 
of full ionization allows quick formation of II2, which then cools the central region before it 
reaches the escape velocity. 
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For this experiment, we initialized ionized fractions as following: ?/hi = 6.4 x 10 ^, a; = 
1.15, ymi = 1, yuei = 6.8 x 10"'^, yndi = 8.9 x 10"^, yueiii = 7x 10~^, yi = for other species. 
Without disturbing the halo density profile - we use the TIS halo model, which is described 
in Section 13.2.11 -, we also assigned a high initial temperature appropriate for photoionized 
gas, T = 2 X lO'^K. These abundance and temperature values roughly mimic the condition 
found in typical H II regions. 

We find that such an initial condition leads to the collapse of the core region, when the 
formation of II2 stimulated by the high initial electron fraction enables H2 cooling. Gas in 
the outskirts evaporates from the halo, however, because pressure forces accelerate the gas 
to escape velocity before it can form H2 and cool. The II2 cooling and adiabatic cooling 



which happen later in this outflowing gas do not rever se the evapora tion 



^ig-ED- 



Our results for this case agree with the outcome of 10 'Shea et al.l (120051 ). This led those 
authors to suggest that the first stars exerted a positive feedback effect on their surroundings, 
triggering a second generation of star formation. A question arises, however, as to whether 
this fully-ionized initial condition of nearby minihaloes is actually achieved b y the first Pop 



Alvarez et al 



III star to form in their neighbourhood. As already mentioned in Section [T], 
(j2006a ) found that the I-front from this Pop III star gets trapped in those minihaloes and 
cannot reach the central region before the sta r dies. In this paper, we will confirm that the 



fully-ionized initial condition of 



O'Shea et al 



(I2OO5I ) is never achieved when one considers 



the coupled radiative and hydrodynamic processes more fully. We will also show that, if any 
protostellar region is to form in the target halo, it does so in the neutral core region which 
the ionizing photons do not penetrate. 



5 MINIMUM HALO MASS FOR COLLAPSE: THE CASE WITHOUT 
RADIATIVE FEEDBACK 

When a minihalo forms as a nonlinear, virialized, gravitationally-bound structure out of 
the linearly perturbed IGM, a change of chemical abundance occurs due to the change of 
gas properties. Most importantly, the hydrogen molecule fraction changes from the IGM 
equilibrium value, yu^ ~ 2 x 10~^, to a new equilibrium value, yu_2 ^ 10""^. Even with such 
a small fraction, H2 can cool gas to T^.^ ^ 100 K, where represents the temperature 
"plateau" that gas in primordial composition can reach by H2 cooling. 

There exists a minimum collapse mass of minihaloes, Mc^min, above which haloes, in the 
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t=2.5 Myr t=10 Myr 




Figure 3. Evolution of a "flash-ionized" TIS halo - i.e. initially fully-ionized at T = 2 X 10"' K of M = 2 X IQ-^Mq. From left 
to right, each panel corresponds to t=0, 2.5, and 7.5 Myrs. Note that initially (t=0) absent H2 is quickly created and cools the 
central region, while the initially hot gas in the outskirts evaporates before it cools. In each panel, from top to bottom, baryon 
density pi,, pressure p, temperature T, velocity v, electron fraction(x; thick)/neutral fraction(yHi; thin), molecule fraction 3/H2 
and the isothermal Mach number Afj are plotted, respectively. Unless specified otherwise, the subsequent figures will follow the 
convention used in this figure. 

absence of external radiation, can form cooling and collapsing cores within the Hubble time 
at a given redshift. The gap between the H2 cooling plateau temperature, Th2, and the 
minihalo virial temperature, Tvir, given by equation (|T9|) is a useful indicator of the success 
or failure of collapse. For instance, at 2; ~ 20, T^ij. ~ 160 K for M = 2.5 x 10^ M©. As 
Tvir — ! 6ven after gas cools to Thj, it cannot collapse fast enough to serve as a site for 
star formation. On the other hand, Tvir ~ 10^ K for M = 4 x 10^ Mq, and the temperature's 
cooling down to Thj ~ 100 K will make the gas gravitationally unstable, which will lead 



Haiman. Thoul. fc Loeb 



to runa way collapse. This argument is supported by the results of 
(Il996bl ). for example, that collapse can occur only in haloes with T^ir ^ 100 K. 

We model the initial minihalo structure by the TIS model as described in Section 13.2.11 
and let it evolve in the absence of radiation, starting from the IGM chemical abundance and 
minihalo virial temperature (Phase I). We determine Mc^min by the criterion 



CcoU 



where tcou is the time at which the central density reac hes nn 



Abel et al. 



(23) 

10^_cm~^ (the density 



2000), and tu is the 



suitable for initiating three-body H2 formation; see e.g. 
Hubble time at a given redshift. 

We find that Mc,min ^ 7 x 10^ M© at ;z = 20 (see Fig. i]). We have plotted the evolution 
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of minihalo centres in the absence of radiation, where each run star ts from Phase I. This i s 
in rough agreement with Mc^min — 1-25 x 10^ M^, the value found bvlMachacek et aLl (|2001l). 



Yoshida et al 



he d iscrepancy is larger with results by iFuUer Sz Couchman] (l2000l ) and 
(12003 ). where they obtain Mc^min — 7x10^ Mq. The biggest contrast exists with iTegmark et al. 



(119971 ). where they find Mc, 



2 X 10^ Mq at z ^ 20, almost 30 times as large as our 



findings. 

We argue that this discrepancy in minimum collapse mass results primarily from how 
well the minihalo structure is resolved. Unless the centre, which gains the highest molecule 
formation rate due to the highest density, is fully resolved, one could be misled by a poor nu- 
merical resolution such that certain low-mass haloes, which can cool and collapse in reality, 
are in hydrost atic equilibrium i n the s imulation. The resolution becomes poorer in the foUow- 



mg sequence: 



Machacek et al. 



(120011 ). which gives the best agreement with our result, used 
an adaptive mesh refinement (AMR) scheme, resolving baryonic mass down to Mf, ~ 5Mq. 
Such high resolution is suitable to resolve even the central part of the smallest minih aloes 
whose total baryonic mass content is roughly 2 — 3 x 10^ Mq. iFuller fc Couchman! (120001 ) and 



Yoshida et al 



(120031), on the other hand, used the smoothed particle hydrodynamics 



Tegmark et all (119971 ) 



(SPH) 



scheme, using SPH particles of mass ~ 40 — 140 x 10^ Mq. Finally, 
used a uniform top-hat model, where there is no radial variation in gas properties such as 
density and temperature, thus the central region is, in effect, completely unresolved. In addi- 



tion, some of the rates used in 



Tegmark et al. 



(119971 ) were not accurate (IFuller &: Couchman 



2000). 



We believe that Mc^min ~ 7 x 10"^ Mq at z = 20 is close to reality, because our 1-D 
spherical setup is based upon the TIS model which is a highly concentrated structure, and the 
resolution of our code is superior to previous calculation^^- It is not our objective, however, to 
settle the exact value of Mc^min. This estimate is based upon our specific criterion described 
in this section, and is subject to change under different criteria. This may also change if one 



adopts a more lea. 
by accretion (see 



istic halo forniation history to account, for instance, for dynamical heating 



Yoshida et al. 



20031 ) . As the haloes we choose are rather conservatively 
divided into successful collapse (for M > 10^ Mq) and failure (for M < 10^ Mq), agreeing 



with AMR simulation result by iMachacek et al. 



2001 



we shall proceed with our choice of 



After this paper was wr itten a new preprint was posted which is consistent with our description here, finding iVfc^min ~ 

fO^ Mq 



I O'Shea fc Norman 



2009) . 
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Figure 4. Evolution of the centre of TIS minihaloes without radiation. Initially, minihaloes have structure described by the 
TIS model (see Phase I depicted by solid lines in Fig. [2]l with the equilibrium value of primordial IGM chemical abundances. 
We define tcoii as the time to reach jih = 10* cm""^, represented by the horizontal dotted lines. 
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parameter space and see how this fate of minihaloes changes as a result of external radiation 
from a Pop III star. 



6 RESULTS: RADIATIVE FEEDBACK ON NEARBY MINIHALOES BY 
AN EXTERNAL POP III STAR 

As described in Section 13.11 we expose target haloes of different mass to the radiation from 
a Pop III star whose spectrum is approximated as a 10^ K blackbody radiation field and 
whose flux is attenuated by the geometrical factor (^-^^ for different values of i^. In this 
section, we summarize the simulation results for both the Phase I (early irradiation) and 
the Phase II (late irradiation) initial conditions. 



6.1 I-front trapping and photo-evaporation 



In all cases, even in the presence of evaporation, we find no evidenc e of penetration of ioni zing 



radiation into the halo core. This is consistent with the resu 
the H II regions of the first Pop III stars and of 



Shapiro et all (120041 ) and 



ts of 



Alvarez et al. 



Ilievetal 



(2006a 



mm for 



) for 



the encounters between intergalactic I-fronts and minihaloes during reionization. There are 
two main reasons for this behaviour. First, the total intervening hydrogen column density is 
initially high enough to trap the I-front outside the core. Second, the lifetime of the source is 
short compared to the evaporation time. If the source lived longer than the eva poration time, 
the I-f ront would eventually have reached the centre of the halo. In that case, IShapiro et c 
(j2004j ) find that the minihalo gas is completely evaporated. In our problem, however, the 
slow evaporation does not allow the I-front to reach the centre within the lifetime of a Pop 
III star. 

The I-front entering the minihaloes propagates as a weak R-type front in the beginning. 
The I-front then makes the transition to the D-type, after reaching the R-critical state. This 



R-critical state is reached when the I-front velocity vi satisfies the following condition: 



Vl = Ct 2 + (c 



1,2 '^1,1 



,2 nO.5 



(24) 



where cj is the isothermal sound speed, cj = yp/p, and subscripts 1 and 2 represent pre-front 
and post-front, respectively. When the I-front propagates into a cold region (T ^ 10^ K), as 
in our problem, this condition is approximately vi ~ 2ci,2 ~ 20kms~^. In all cases, we find 
that this transition occurs in times less than the lifetime of the source star, 2.5 Myrs. After 
reaching the R-critical state, gas in front of the I-front forms a shock, which then detaches 
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Table 1. Ionized mass fraction of baryons for different mass target fialoes (columns) at different distances from a 120 Mq Pop 
III star (rows; fluxes in units of 10^" s~^ kpc~^, Fq, in square brackets). The ratio shown here is the mass ionized during the 
hfetime of the star to the total baryon mass. 



Total Halo Mass in 10^ Mq units 
(Halo Baryon Mass in 10^ Mq units) 



D (pc) [Fo] 


0.25 
(0.043) 


0.5 ■ 10* 

(0.086) 


1 

(0.17) 


2 

(0.34) 


4 

(0.69) 


8 

(1.371) 


180 [46.3] 


0.95 


0.92 


0.88 


0.84 


0.82 


0.79 


360 [11.6] 


0.85 


0.81 


0.77 


0.74 


0.70 


0.67 


540 [5.14] 


0.78 


0.74 


0.70 


0.66 


0.62 


0.59 


1000 [1.5] 


0.66 


0.60 


0.55 


0.50 


0.47 


0.43 



from the slowed I-front. As an example, we plot in Fig. [5] the profiles of Phase I, 4 x 10^ Mq 
halo at t = tpt-crit under different fluxes. 

All of the post-front (ionized) gas, initially undisturbed, eventually evaporates away, 
accelerated outward by a large pressure gradient. As the line-of-sight is cleared by this 
evaporation, ionizing radiation penetrates deeper, until the source turns off. See Figs \5\ E] 
and [7| for the evolution of the I-front. 

This result invalidates the initial condition adopted bylj 



O'Shea et al. 



(120051 ) and lMesinger et al. 



(120061 ) which led them to find that H2 formed in the core region after it was ionized and then 



cooled while recombining, once the source turned off. As we show, the core remains neutra l 



O'Shea et al. 



(120051) 



before and after the source is turned off, so the mechanism explored by 
does not work. This neutral core, therefore, must find a different way to cool and collapse if 
star formation is to happen in the target minihalo. 

What happens to the initially ionized gas after the star turns off? This gas recombines 
as it cools radiatively and by adiabatic expansion, even forming H2 molecules. We find that 
this cooling cannot reverse the evaporation, however. Gas is simply carried away with the 
initial momentum given to it when it was in an ionized state. In Tabled], we list the fraction 
of the baryonic halo mass which is ionized during the lifetime of the star. This mass serves 
as a crude estimate of the mass lost from these haloes by evaporation. We found no major 
difference between Phase I and Phase II in this matter, so we provide only one table. 



6.2 Formation of H2 precursor shell in Front of the I-Front 

We find that a thin shell of H2 is formed just ahead of the I-front, with peak abundance 
?/H2 ~ 10~^. It happens mainly because the increased electron fraction across the I-front 
promotes the formation of H2. More precisely, the gas ahead of the I-front is ionized to the 
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Figure 5. Radial profiles of 4 X 10^ Mq fialo at t = tR_crit in Phase I, for different fluxes (distances). R-critical condition is met 
by tlie condition vi si 2cs i 2 ~ 20kms~^, when the I-front makes a transition from R-type to D-type. Different distances to the 
source are represented by different line types: D = ISOpc {Fq = 46.3; black, solid), D = 360pc {Fq = 11.6; red, short-dashed), 
D = 540pc [Fq = 5.14; green, dotted), D = lOOOpc (Fq = 1.5; blue, long-dashed). In the fifth panel from top, electron fraction 
X and neutral hydrogen fraction j/jj (thin) are plotted. As is shown in the density (pt,) plot, gas just starts to respond to the 
I-front hydrodynamically, as the initial R-type, supersonic I-front slows down to reach the R-critical phase. Also note that the 
shorter the distance, the deeper the I-front is (i.e. smaller r) at i = ta_crit- 

extent that the electron abundance is large enough to form H2, but at the same time too low 
to drive significant coUisional dissociation of H2. The width of this H2 shell and the amount 
of H2 in this region is determined by the hardness of the energy spectrum of the source: the 
width of the I-front is of the order of the mean free path of the ionizing photons. Pop III 
stars, in general, produce a large number of hard photons due to their high temperature, 
which can penetrate deeper into the neutral region than soft photons. 

This precursor H2 shell feature is evident in Figs [5], [6l and [3 We show the detailed 
structure of these H2 shells in Fig. [8], where we plot the radial profile of the abundance of 
different species for the case of M = 4x 10^ Mq, Phase I, D = 540 pc (Fn = 5.14) at t = 0. 5 1^. 



Ricotti. Gnedin. fc Shulll (120011) for 



Ricotti et al 



2001) 



Susa fc Umemura 



We note the similarity between our results and those of 
an I-front in a uniform, static IGM at the mean density (see Fig. 3 in 
which also show a precursor H2 shell. A similar effect was reported by 
mom , as well. 

What is the importance of this II2 shell in protecting the central region of haloes from 
dissociating radiation? The molecular column density obtained by this H2 shell sometimes 
reaches ~ 10^^ cm~^, which provides an appreciable amount of self-shielding. The self- 
shielding due to the H2 shell, however, is not the major factor that determines whether 
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es of 4 X 10^ Mq halo at t = 0.5t» in Phase I, for different fluxes (distanees). Same line types are used 




Figure 7. Radial profiles of 4 X 10^ Mq halo at t = t* in Phase I, for different fluxes (distances). Same line types are used as 
in Fig. |5] 



or not the H2 in the core region is protected. A more important factor is which evolutionary 
phase the target halo is in when it is irradiated. Roughly speaking, when a target halo is 
irradiated early in its evolution (Phase I), the precursor H2 shell dominates the total H2 
column density available to shield the central region, but this shielding is not sufficient to 
prevent photodissociation there anyway. On the other hand, if the halo is irradiated later in 
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r (pc) 



Figure 8. Radial profile of abundance of primordial gas species at t = 0.6 1», for a halo of M = 2 X 10^ Mq in Phase I, 
illuminated at D = 540 pc {Fq = 5.14). Labels are self-explanatory; electron abundance x is represented by the dotted curve, 
which closely follows the H II abundance. Th e flux is coming from the right hand side, so this figure can be compared to the 
mirror image of Fig. 3 in lRicotti et al. 1 1I2OOII) for an I-front in a static, uniform IGM at the mean density. 

its evolution (Phase II), the H2 column density of the shell is only a small part of the total 
H2 column density, so shielding is successful independent of the precursor shell. We describe 
this in more detail as follows. 

In order to understand quantitatively the importance of the H2 shell in protecting the 
central H2 fraction, we have performed simulations with a source SED that is identical to 
the Pop III SED below 13.6 eV, but zero above 13.6 eV. As the radiation is now incapable of 
ionizing the halo gas, the H2 shell formation by partial ionization will not occur. This enables 
us to compare our results where the H2 shell is present to those cases without an H2 shell. 
We describe a specific case of M = 2 x 10^ Mq as an illustration. Roughly speaking, the H2 
shell which forms only in the presence of ionizing radiation compensates for the amount by 
which the initial molecular column density, A^H2; is reduced when molecules in the ionized 
region are destroyed by collisional dissociation. The nett column density in the case where 
the H2 shell is present even exceeds that in the case without the H2 shell (Figs [9] and [TOl) . 
The nett effect is the increase of the self-shielding. Such an increase of the self-shielding, 
however, is not too dramatic. In the case of M = 2 x 10^ Mq with Phase I initial conditions, 
?/H2 ~ 10"^ '^ at the centre, about an order of magnitude higher than the central 1/H2 of 
the case without ionizing photons (Fig. [9]). This molecule fraction is still too low, however, 
to cool the gas. On the other hand, in the case of M = 2 x 10^ Mq with Phase II initial 
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Figure 9. Radial profiles of H2 fraction (top) and the H2 column density (bottom) for M = 2 X 10^ Mq, Phase I initial 
conditions, at t = (dotted) and at t = t» (solid). The source is at a distance D = 540 pc (_Fo = 5.14). Also plotted are 
those for a radiation composed only of dissociating photons (dashed) at t = t» . Even though H2 shell provides self-shielding by 
contributing A'jjj ~ 10^® cm^, the molecule fraction yn^ is held at ^10~^'^ due to strong dissociating radiation. 



conditions, y^^ ~ 10~'^'^ at the centre throughout the hfetime of the Pop III source, whether 
or not the H2 shell is formed. The depth (radius) of penetration of dissociating photons 
differs by a factor of 2 if the shell is included, but the central H2 is still protected because of 
the high H2 column density apart from the precursor shell (Fig. [TOl) . The major factor that 
determines the fate of the central H2 fraction is instead the evolutionary phase of a target 
halo when it is irradiated. The short lifetime of a Pop III star plays an important role of 
either reconstituting or protecting molecules in the core, depending upon the evolutionary 
phase of the halo, as will be described in Section 16.41 

Note that in all cases, we use equation (I13p . the shielding function for thermally-broadened 
lines with T = 10^ K. This is justified by the fact that the H2 shell moves inward with 
V ^ 2 — 5kms~^ and the shell achieves T = Tgh ~ 10'^ — 5 x 10'^ K, where Tgh denotes the 
temperature of the shell. If we take this peculiar velocity as sound speed, v ^2 — 5kms^^ 
corresponds to T = Tp = v'^firriYi/k = 6 x 10^ — 3.7 x 10'^ K, where the subscript p denotes 
the peculiar velocity. A crude way to imitate both effects by thermal broadening is to use 
the sum of these two temperatures (Tsh and Tp). We take the most conservative stand - the 
least self-shielding effect - in order not to overestimate the self-shielding, and use T = 10^ K 
as the temperature responsible for the nett thermal broadening of the molecular LW bands. 
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Figure 10. Radial profiles of H2 fraction (top) and the H2 column density (bottom) for M = 2 X 10^ Mq, Phase II initial 
conditions, at t = (dotted) and at t = 0.5t, (solid). The source is at a distance D = 540 pc {Fq = 5.14). Also plotted are 
those for a radiation composed only of dissociating photons (dashed) at i = t» . Contrary to Fig. [O] the core H2 fraction is well 



protected to a high level, yn 



IQ- 



The role of H2 shell is, however, not dramatic, because even without this shell, the 



core is protected from dissociating radiation (dashed). It rather offsets the loss to the molecular column density from coUisional 
dissociation in the ionized region. 



6.3 Formation of shock and Evolution of core 

After the I-front decelerates as it enters the target halo, transforming from R-type to D-type, 
a shock front forms to lead the D-type front. The neutral gas in the core is strongly affected 
by this shock front as it propagates. This shock plays an important role in providing both 
positive and negative feedback effects. By identifying successive evolutionary stages of the 
shock, we now describe how the core responds to the shock and evolves accordingly. 



6.3.1 Stage I: Formation and acceleration of Shock 

A shock starts to form as the I-front, initially moving supersonically as an R-type, slows 
down and turns into a D-type. The pre-front gas - neutral gas ahead of the I-front - can 
respond to the I-front before it is swept by the I-front, because the D-type front moves 
subsonically into the neutral gas. It is easier to understand the formation of the shock by 
using the I-front jump conditions: the pre-front gas speed in the rest frame of the I-front, Vi, 
derived from the I-front jump conditions, should satisfy either f i > fj? = ci^2 + (cj 2 ~ Cj ^)°'^, 
or vi < vd = ci,2 — (ci2 — Cj i)"'^, where ci^i and ci,2 are the isothermal sound speeds of the 
pre-front and post-front gas, respectively, vji and vd have a gap of 2(cj2 — Cj which 
is nonzero in general. As the I-front slows down and Vi starts to cross vr, Vi encounters a 
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value which is not allowed mathematically. This paradox is resolved, however, because the 
pre-front gas now "prepares" a new hydrodynamic condition by forming a shock. The shock 
wave increases pi and thereby reduces vi and increases vd, making it possible to satisfy the 
D-type condition, vi < Vfj. 

This shock-front then propagates inward, separating from the 1-front, due to the discrep- 
ancy between the speed of the shock-front and the speed of the I- front. As the shock-front 
enters the fiat-density core, the shock front starts to accelerate, leaving behind the post- 
shock gas with ever increasing temperature (e.g. see time steps 4 and 5 in Fig. [TTl where 
the post-shock temperature increases as the radius r decreases). 

As the shock boosts the density and temperature in the neutral, post-shock gas, the 
H2 formation rate there increases, boosting the H2 column density even further. We can 
understand the evolution of 1/H2 the presence of this shock quantitatively by using its 
equilibrium value, ?/H2,eq- The increase of density and temperature due to this shock pro- 
motes H2 formation, as follows. When there is no significant destruction mechanism, 
the dominant H2 formation mechanism is through H~ (equation [T]), and the H2 formation 
rate becomes equivalent to the H~ formation rate. Photo-dissociation dominates over coUi- 
sional dissociation in destroying H2, which occurs when x < 4 x lO^^Tj^^^ and nn > 0.045 x 



coefficient (jde Jong)ll972l ) 



(Flw/10 erg s ^ cm ^ H z ^) (e.g. iGlover fc Brand! 1200 ll ). Using the H formation rate 



= 10"^^ Tk cm^ s"\ (25) 
and the photo-dissociation rate coefficient k^^ given by equation ( fT2l) . we obtain 

where we have used the fact that one can scale Flw by Fq according to the following: 

Flw ~ 3.25 X 10-2^ erg s"^ cm'^ Hz"^ Fq, (27) 

if one adopts a black-body spectrum with T = 10^ K. As seen in equation (l26l) . both the high 
temperature (~ 1000 — 5000 K) and increased density (x4 in the case of strong shock) of 
the post-shock gas contributes to boosting the H2 fraction. As yn^ oc -Fshieid; niolecular self- 
shielding also plays an important role in determining ?/h2- If the shock boosts the formation 
rate of H2 and yn^ increases, so will A^H2) with it the shielding. These two effects, 
therefore, amplify each other. 
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There is an additional mechanism to create molecules: the shock-induced molecule forma- 
tion (SIMF). The acceleration of the shock-front accompanied by an increasing post-shock 
temperature, leads to a partial ionization of the post-shock gas in many cases, when the 
right condition (T > 10^ K) is met to trigger collisional ionization - see, for example, step 5 
in Fig. [TTJ the centre is shock-heated above 10^ K, with a boost in x. The electron fraction 
X now reaches ~ 10~^ — 10"^, which promotes further H2 formation. This mechanism is 
indeed identical to the formation mechanisrn in a gas that has be en shock-heated to 



temperatures above 10^ K ( IShapiro &: Kanglll987l : iKang fc Shapirol 19921 ). When a gas cools 



radiatively from a temperature well above 10^ K, it cools faster than it recombines. As a 
result, the recombination is out of equilibrium, and an enhanced electron fraction exists at 
temperatures even below 10^ K compared to the equilibrium value. This electron fraction 
triggers the formation of H2 through the gas-phase reactions (equations [1] and [2]). 

SIMF does not always occur, however. The shock-front can accelerate when the pre-shock 
density remains almost constant (e.g. Fig. [TTl) . If the density increases faster than the shock 
propagates, on the other hand, the shock-front will encounter an ever increasing density 
"hill" and it will never accelerate to generate post-shock temperature above 10^ K (e.g. Fig. 
T2I) . The dependence of SIMF on the halo mass, source flux, and the initial phase will be 
described in Section [6^ 



6.3.2 Stage II: Cooling and Compression of Core 

As the shock-front approaches the centre of the halo, the post-shock gas there becomes more 
concentrated and denser than the pre-shock gas. This shock-induced compression leads to a 
very fast molecular cooling in the core and further compression in almost a runaway fashion, 
as follows. 

Molecular cooling occurs very rapidly at a high density and temperature condition. As- 
suming that the pre-shock gas of the halo core remains unchanged before the shock-front ar- 
rives - as is usually the case in Phase I - and the shock is strong, the post-shock density of the 
core becomes 4 times higher than that of the pre-shock, namely rini ~ 4x30 cm~^ = 120 cm~^ 
in a TIS halo core at 2; = 20. At the same time, post-shock temperature can be as high as 
10^ K. The molecular cooling time, tcooi,H2 = T/{dT/dt), is 
kT 



''cool,H2 



X/i(7 - 1)?/h2^hiAh2 
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where X = 0.75 is the hydrogen mass fraction, and is the molecular cooling rate. For a 
gas with nui = 120 cm~^ and T = 10^ K, Arj ~ 3.4 x 10~^^ erg cm~^ s~^, and thus 

tcooi,H2 ~ 1.8x 10=^yr (^-^^ \ (29) 

With such a rapid cooling, the isothermal shock jump condition (T2 = Ti) is a good approx- 
imation, and the post-shock density becomes even higher than that of the adiabatic strong 
shock, because Pb,2/Pbi ~ Mf^ now. Such a strong compression of the core is observed very 
frequently in our parameter space of different halo masses and source fluxes. For example, 
Fig. [I3] shows how the centre of a halo with M = 2 x 10^ Mq evolves in response to the 
shock. As the shock hits the centre, density increases by many orders of magnitude. 

Does this compression eventually lead to the core collapse? As the shock carries the 
kinetic energy as well as the thermal energy, the shock will bounce off the centre after it hits 
the centre. In the following section, we describe this final stage of the shock propagation and 
show how it will affect the core collapse. 



6.3.3 Stage III: Bounce of Shock and Collapse of Core 

After the shock hits the centre, the shock wave will be reflected and propagate outward. In 
our ID calculation, this reflection will mimic the transmission of the shock wave through the 
centre. This bouncing shock will try to disrupt the gas. The core that is undergoing cooling 
and compression due to the positive feedback effects mentioned so far will be affected by 
this negative feedback effect, as well. 

The final fate of the core depends on how well the core endures such a disruption. As 
the shock bounces off the centre, density starts to decrease. If this bounce is weak, the core 
quickly reassembles, cools, and finally collapses. If this bounce is strong, the core will take 
a longer time to collapse and, in some cases, the core will never collapse within the Hubble 
time. Haloes of smaller mass seem to be more susceptible to this shock-bounce than those 
of larger mass (see Figs and [T5l for comparison). 

If the core finally takes the collapse route, the central hydrogen number density increases 
to ~ lO^cm"^, at which point the ro-vibrational levels of H2 are populated at t heir equilib- 



rium values and the molecular cooling time becomes independent of density (e.g. 



Abel et al. 



2OO2I ). Since then, adiabatic heating dominates over the molecular cooling, and the temper- 
ature increases as collapse proceeds. Finally, when nui reaches ~ lO^cm"^, the three-body 
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hydrogen reaction ensues and converts most hydrogen atoms into the molecules, which will 
undergo a further collapse and form a proto-star. 

6.4 Feedback of Pop III starlight on Nearby Minihaloes: parameter 
dependence of core collapse 

We now summarize the outcome of our full parameter study of radiative feedback effects 
of Pop III starlight on nearby minihaloes. As we have described in the previous section, 
positive and negative feedback effects of the shock compete and produce a nett effect which 
can be either 1) an expedited collapse, 2) delayed collapse, 3) neutral (unaffected) collapse, 
or 4) a disruption. 

Overall, the radiative feedback effect of a Pop III star is not as destructive as naively 
expected. Minihaloes with M> [1 — 2] x 10^ Mq, which can cool and collapse without radia- 
tion, are still able to form cooling and collapsing clouds at their centre even in the presence 
of Pop III starlight. The quantitative results are summarized in Tables IH [3] and Fig. [T6l 

The relatively short lifetime of a Pop III star, compared to the recombination timescale 
in the core, is a key to understanding this behaviour. One of the necessary conditions for 
the core collapse is that H2 molecular cooling should occur in the core. As this requires a 
sufficient molecular fraction, namely > 10"'^, it is crucial to understand how molecules 
are created at such a level. In Phase I (low y-^^ and high x), radiation can easily dissociate 
H2 while the source is on, but after the source dies, the high electron fraction stimulates H2 
formation. This is possible because the recombination time in the TIS core is longer than 
the lifetime of the source Pop III star. On the contrary, in Phase II (high J/H2 ^^^id low x), 
B.2 is more easily protected against the dissociating radiation because the higher H2 column 
density provides self-shielding and compression increases the formation rate. Because the 
source irradiates these haloes for a short period of time, the dissociation front does not 
reach the centre, and its high molecule fraction is preserved throughout the Pop III stellar 
lifetime. 

6.4.1 Phase I 

When haloes start their evolution from Phase I - IGM chemical abundance and the TIS 
structure -, other than the change of collapse times, there is no reversal of collapse. In other 
words, haloes that were destined to cool and collapse would do so even when exposed to the 
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M=2e5 Mg; D=540 pc M=2e5 M^; D=540 pc 




Ig r (pc) Ig r (pc) 

Figure 11. Success of shock-induced molecule formation (SIMF): evolution of 2 x 10^ -Mq halo irradiated early (Phase I) by a 
Pop III star at D = 540pc (_Fo = 5.14). Radial profiles of fluid parameters — baryon gas density (pt), pressure (p), temperature 
(T), radial velocity (v), electron fraction (x; thick line), neutral fraction {{/hi! thin line), and molecule fraction {yn^) ^ ^re 
labelled by different time frames as following: 0-t = 0;l-t = tn-crit = 0.2 1,; 2 - t = 0.5 1,; 3-t = t,;4-t = 1.5 1,, 5 - 
* = *shock bounce = 1.611 1,; 6 - t = tshock bounce + £ = 1.617 1,; 7-t = 2t, ;8-t = tcoll = 2.6 1*. These time frames arc shifted 
along the y-axis for clarity, with equal displacements as following: Algp = 1; Algp = 1; AlgT = 2; Av = 5 km/s; Alga; = 10; 
AlgJ/Hi = 10; Alg2/H2 = 5- E)otted lines represent the initial central density, T = 10* K, v = Okm/s, and = 10~* in the 
Pl,, T, V, and plot, respectively, tshock- bounce is the time when the shock front reaches the centre. Note that at this moment 
the shock-front accelerates to heat the gas up to T > lOOOOK at the centre. This temperature is high enough to cause coUisional 
ionization, which leads to rapid formation of H2 and cooling at the centre afterwards (t ^tshock-bounce +£)• Thus the thermal 
energy delivered is dissipated very easily, and the core collapses in a runaway fashion. We show here the fast evolution of H2 
around the time of shock-bouncing, using e f» 1.5 x 10* years. 
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M=2e5 Mg; D=540 pc M=2e5 M^; D=540 pc 




Ig r (pc) Ig r (pc) 

Figure 12. Failure of slioclt-induced molecule formation (SIMF): evolution of 2 X IO^Mq halo irradiated late (Phase II) by 
a Pop III star at D = 540pc (_Fo = 5.14). Contrary to the case with the Phase I initial conditions with the same mass and 
flux (Fig. [TTI l. SIMF does not occur in this case. Gas profiles are labelled by different time frames as following: - t = 0; 1 - 
t = tR_crit = 0.2U; 2 - t = 0.5 1*; 3 - i = t,; 4 - t = 1.5t,; 5 - t = tshock bounce = 1.67t*; 6 - t = 2t*; 7 - i = Wl- These 
time frames are shifted in the same way as in Fig. 1111 Dotted lines have the same meaning as those in Fig. 1111 Note that even 
at t = tshock— bounce 1 the shock-front velocity is not high enough to heat the gas up to T> lOOOOK at the centre. The SIMF, 
therefore, does not occur. The thermal energy delivered, however, is dissipated anyway by radiative cooling, because the core 
is well protected from the dissociating radiation and the high H2 fraction is maintained throughout the evolution. The core 
collapses in a runaway fashion afterwards. 



first Pop III star in the neighbourhood. Minihaloes with M > 10^ Mq are able to collapse 
without radiation, while those with M < 10^ Mq are not. In the presence of radiation, haloes 
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Figure 13. Evolution of the central part of a Af = 2 X 10^ Mq halo, for Phase I (left) and Phase II (right) at D = 180 pc 
(top; Fq = 46.3) and D = 540 pc (bottom; Fg = 5.14). In Phase I, expediting core collapse is observed for both distances. In 
Phase II, mixed results occur: delayed collapse for D = 180 pc while expedited collapse for D = 540 pc. Another notable feature 
is the shock-ionization (electron) molecule formation in all cases except for the case of D = 540 pc. Phase II. The increase of 
molecule fraction in the latter case is due to the increase of temperature and density due to shock compression, wrhile in other 
cases, shock-induced electron formation promotes further molecule formation. 



with M > 10^ Mq are still able to collapse, while those with M < 10^ Mq are still unable to 
do so, even with the help of shock-induced molecule formation (Fig. [T6j Table [2]). 

The core collapse in Phase I occurs mostly as an expedited collapse (Table [2]). The shock 
plays a major role in driving such an expedited collapse: the H2 fraction becomes boosted 
by the higher density and high temperature delivered by the shock. Whether or not SIMF 
has occurred, such a boost in is sufficient to expedite the core collapse. 

There is one delayed collapse case at the low mass and the high flux end. For M = 10^ Mq 
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Figure 14. Evolution of the central part of a M = 1 X 10''' Mq halo, for Phase I (left) and Phase II (right) at D = 180 pc 
(top; _Fo = 46.3) and D = -540 pc (bottom; Fq = 5.14). In Phase I, high flux results in the delayed collapse (D = 180 pc), 
while low flux results in the expedited collapse {D = 540 pc). In Phase II, core collapse is completely halted at any flux. The 
shock-induced molecule formation occurs in all cases, but the negative feedback is stronger than the case of higher masses. 



at Fo = 46.3, the boosted molecule formation is not sufficient to bring the core to an 
immediate collapse. As the shock bounces, the momentum carries gas away from the centre 
until it cools and recoUapses. 

The unchanged collapses occur at the high mass and the low flux end. For M — 8x 10^ Mq 
at Fq = [1.5, 5.14], the shock propagates into the already collapsing core. The shock energy 
dehvered in these cases is not significant enough to change the course of collapse. 
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Figure 15. Evolution of the central part of a M = 8 X 10'' Mq halo, the highest mass end of our parameter space, for Phase I 
(left) and Phase II (right) at D = 180 pc (top; Fg = 46.3) and D = 540 pc (bottom; Fq = 5.14). An expedited collapse occurs 
for high flux {D = 180 pc). Otherwise, collapse is unchanged in time. 



6.4.2 Phase II 

The overall effect of radiation from a Pop III star on neighbouring minihaloes in Phase II 
is similar to the effect on the minihaloes in Phase I: haloes that were destined to cool and 
collapse would do so even when exposed to the first Pop III star in the neighbourhood. A 
slight shift of the trend exists, however, in Phase II (Fig. [161 Table [3]). When haloes start 
their evolution from Phase II, those with M > 10^ Mq are able to collapse without radiation, 
while those with M<2 x 10^ Mq are not. The collapse in Phase II is reversed (halted) for 
the low mass end: for M = 10^ Mq, the shock disrupts the core and it never recollapses. 
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Figure 16. Feedback of Pop III starlight on Nearby Minihaloes: To collapse or not to collapse? Outcome depends on minihalo 
mass and stellar flux Fa = /-'/(lO s^^ kpc^^) (oc distance"^) as plotted and on the timing of the feedback. The flux (distance) 
is 46.3 (180 pc), 11.6 (360 pc), 5.11 (540 pc) and 1.5 (1000 pc), from top to bottom in each panel, for a 120 Mq Pop 111 star 
located at each distance. The two panels correspond to different initial conditions when starlight arrives: (left) (Early = Phase 1), 
virializcd halo in hydrostatic equilibrium with IGM primordial chemical abundances and (right) (Late = Phase 11), halo evolved 
chemically and hydrodynamically without radiation until x = 10^'' at centre. The outcome of the radiative feedback is marked 
by a circle (collapse) or a triangle (no collapse), as well as the logarithmic collapse time (size of circle). Solid dots represent 
those cases in which shock-induced molecule formation (SIMF) occurs. Compared to the "no radiation" cases on the bottom of 
each panel, feedback malses collapse either (1) expedited (smaller circle), (2) delayed (larger circle), (3) unchanged (same sized 
circle), (4) reversed/failed (triangle), or (5) unchanged/no collapse (triangle). Expedited or unchanged collapse occurs widely 
for 10^ < M/Mq<8 X 10^ in the left panel (Phase I), with exception of M = 10^ Mq and Fq = 46.3(1? = 180 pc). For right 
panel (Phase II), reversed cases occur for M = 10^ Mq, delayed collapse for M = 2 x 10^ Mq and Fo = 46.3(13 = 180 pc), and 
expedited or unchanged collapse for the rest. 
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Table 2. Collapse times of Phase I for different target haloes (columns) at different locations (rows). Each element represents 
the ratio tcoii,R/*coU,NR> where tcoll.R is the collapse time (time the halo core takes to reach n^rit = 10®cm~^) under radiation, 
and tcoll,NR the collapse time without radiation. tcoll.NR is denoted by values in parentheses. Dot represents the case where 
the core collapse never orcurs (hiring tho Hnl)})Ie tiiiio at z = 20. or 186 million A'cais. 



Total Halo Mass in lO"' Mq units 
(Collapse Time without Radiation in Myrs units) 



D (pc) [Fo] 0.25 


0.5 1 


2 




4 


8 


(•) 


(•) (88.82) 


(31.02) 




(14.61) 


(8.66) 


180 pc [46.3] 


1.455 


7.288 ■ 10" 


2 


1.838 ■ 10-1 


4.712 ■ 10-1 


360 pc [11.6] 


1.935-10-1 


1.308 ■ 10- 


1 


3.597- 10-1 


8.177 ■ 10-1 


540 pc [5.14] 


3.427 ■ 10-1 


2.093 • 10" 


1 


4.919 - 10-1 


1.000 


1000 pc [1.5] 


9.497-10-1 


4.525 ■ 10" 


1 


7.144- 10-1 


1.241 



Table 3. Collapse times of Phase II for different target haloes (columns) at different locations (rows). Each clement represents 
the ratio tcoll,R/tcoll,NR) where tcoll,R is the collapse time (time the halo core takes to reach n^rit = 10** cm"'^) under radiation, 
and tcoll.NR the collapse time without radiation. icoll.NR is denoted by values in parentheses. Dot represents the case where 
the core collapse never occurs during the Hubble time a,t z = 20, or 186 million years. 



Total Halo Mass in 10^ Af© units 
(Collapse Time without Radiation in Myrs units) 



D (pc) [Fo] 0.25 0.5 


1 2 


4 


8 


(■) (■) 


(65.66) (14.49) 


(4.23) 


(1.65) 


180 pc [46.3] 


4.269 


7.151 ■ 10-1 


9.541 ■ 10-1 


360 pc [11.6] 


4.997 ■ 10-1 


1.155 


1.002 


540 pc [5.14] 


6.740-10-1 


9.794 - 10-1 


9.964-10-1 


1000 pc [1.5] 


5.794-10-1 


9.926 - 10-1 


9.994-10-1 



SIMF occurs at Fq > 1.5 for M — 10^ M©, but this does not prevent such a destructive 
process from happening. 

As haloes start their evolution from Phase II, in which the halo cores are already cooling 
and collapsing, the neutral (unaffected) collapse cases occur more frequently than in Phase 
I. At high and intermediate masses, the collapse time hardly changes from the case without 
radiation. Haloes with M — 8 x 10^ Mq collapse before the source dies, as they do with- 
out radiation, simply because the shock wave does not affect the core. In this case, shock 
propagates into the centre after collapse has advanced significantly. 

There is one delayed collapse case: compared to the delayed collapse in Phase I, which 
occurred at low mass/high flux end (M = 10^ M© at Fq — 46.3), this now occurs at an 
intermediate mass/high flux end (M — 2x 10^ Mq at Fq = 46.3). Otherwise, for intermediate 
mass, collapse is either neutral or expedited. 
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6.5 The structure of haloes at the moment of collapse 



The structure of halo at collapse determines how a protostar evolves into a star and how 
the starlight will later propagate through the host halo. We first show how halo profiles at 
collapse vary for different mass without radiation. We then describe how halo structure is 
affected by the Pop 111 starlight. 

We note that halo structure shows a strong dependence on the halo mass. For radius 
r > 10^^ pc, density profiles of haloes without radiation are well fit by a power law, p oc r~^. 
The value of w, however, is dependent upon the mass of the halo. We find that w = 2.5, 2.4, 
2.3, and 2.2 for haloes of mass M = 10^ M = 2 x 10^ M = 4 x 10^ and M = 8 x 10^ M©, 
respectively. In all cases, the temperature is somewhat fiat with T ~ 10^'^ — 10^ K. The 
temperature at r ^ 10~^pc, where p ^ 3 x 10~^^gcm~'^ (or rin ^ lO^cm"^), is about 800 K 
in all cases. The universality of these core properties seems to originate from the fact that 
the dominant process, H2 cooling, causes loss memory of the initial condition (e.g. different 
virial temperatures for different virial masses). The outer part of these haloes, however, still 
retain the memory virial equilibrium because radiative cooling is negligible. Overall, as mass 
decreases, density slope increases (see Fig. [T7I) . 

The radiative feedback effect of the starlight on final halo profiles is found to be negligible 
in most cases. The region that has been photo-ionized during the stellar lifetime is obviously 
strongly affected. The neutral region, however, is almost indistinguishable from the case 
without radiation in most cases. The variance of temperature profile exists only at the 
low-mass end, M = 10^ M©, or the high-fiux end, Fq = 46.3 {D = 180 pc). Such variance 
completely disappears at the high-mass end, M = 8 x 10^ Mq, because collapse is mostly 
unaffected (Fig. [T7|l . 

This result indicates that the mass of secondary Pop III stars would be almost identical to 
that of the Pop 111 stars which form without radiative feedback effect. A more fundamental 
variance may exist , howe ver, due to the environmental variance of star forming regions: 



O'Shea fc Norman! (120061 ) show that temperature variance of different regions result in the 
variance of protostellar masses, due to the corresponding variance of mass infall rate. As our 
simulation does not advance beyond nn = 10^ cm~'^, where three-body collision can produce 
copious amount of H2 molecules and change the adiabatic index of the gas, we are unable 
to quantify the final mass of the protostar at this stage. 
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Figure 17. Halo profiles at the onset of collapse. The top left panel shows profiles of different mass haloes that collapse without 
radiation, with mass M = 10^ Mq (solid), M = 2 X 10^ Mq (dotted), Af = 4 X 10^ Mq (short-dashed), and M = 8 X 10^ M© 
(long-dashed). Other panels show profiles of haloes of different masses and phases. In each of these panels (except for the 
top-left panel), no radiation (black; solid), D = 180 pc {Fg = 46.3; red; dotted), D = 360 pc (_Fo = 11.6; green; short-dashed), 
D = 540 pc (_Fo = 5.14; blue; long-dashed), and D = 1000 pc (_Fo = 1.5; cyan; dot-dashed) cases are plotted. Note that even 
though the region ionized during the stellar lifetime is heavily affected, the final structure of the neutral core remains hardly 
changed at the onset of collapse in most cases. Some variation is seen at the low-mass end, M = 10^ Mq, or the high-flux end, 
D = 180 pc (Fo = 46.3). 
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6.6 Feedback of Pop III Starlight on Merging Haloes and Subclumps 

While we were preparing this manuscript, two preprints were posted describing simulations 
of the radiative feedback of the first Pop 111 star on dense gas clumps even closer to the 
star than the external minihaloes we have considered so far, for the case of subclumps 



( ISusa fc Umemurall2006l ) and the case o: 



a second miniha lo undergoing a major merger with 



the minihalo that hosts the first star ( lAbel et al.l l2006l ). The centre of the target halo or 



clumps in this case is well within the virial radius of the halo which hosts th e first star, 



and, t hese authors find that secondary star formation occurs in these subhaloes. 



Abel et al. 



( I2OO6I ). for instance, report that the first star forms inside a minihalo of mass M = 4 x 10 Mq 
as it merges with a second minihalo of mass M = 5.5 x 10^ Mq (the target halo). The centre 
of this target halo is at a distance of only 50 parsecs from the first star. Cooling and collapse 
leading to the formation of a protostar is found to occur inside the target halo about 6 Myrs 
after the first star has died. 

We ask the same question that whether or not a halo would collapse to form a secondary 
Pop 111 star if a nearby Pop 111 star irradiates the halo at a distance of 50 pc. Note that the 
target halo we consider now would collapse anyway if there were no radiation, in ~ 1 1 Myrs 
for Phase 1 and ~ 3 Myrs for Phase 11 (see Table HI). This problem requires us to extend 
our parameter space beyond what has been considered so far, because of the short distance 
(high flux) between the source and the target. 

We have attempted to reproduce the result of 



Abel et al. 



(120061 ) using our code for a 



target halo of mass M = 5.5 x 10^ and D = 50 pc, corresponding to the ionizing flux 
Fo = 600. Note that the LW band flux is very high: Flw ~ 2000 x 10"^^ ergs"^ cm^^ Hz^^ 
(equation [27|) . As D is smaller than the virial radius of the target halo, we truncated the halo 
proflle at 50 pc. To be consistent with our previous calculations, we neglect the geometrical 
variation of the flux with position inside the target halo. 

Surprisingly enough, contrary to the outcome of L 



Abel et al. 



(120061 ). we flnd that collapse 



is expedited, occurring within the lifetime of the first star, for both Phase 1 and Phase II 
initial conditions. The main mechanism was SIMF: initially, H2 is completely wiped out by 
a strong dissociating radiation, but as the SIMF occurs, newly created molecules l ead t o 
cooling and collapsing. This result is in disagreement with the result of Abel et al. ( 20061 ). 
which shows that the second star forms after the star has died. 

This puzzling result shows the importance of H2 self-shielding. 
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formed an optically-thin calculation for Lyman- Werner bands, neglecting the H2 self-shielding 
while our calculation took the self-shielding into account. In order to mimic their calculation 
more consistently, we artificially performed an optically-thin calculation for Lyman- Werner 
bands. We found that, if the target halo is irradiated without H2 self-shielding, the core 
collapse is delayed and occurs after the star dies both in Phase I and Phase II. In our simu- 
lations without H2 self-shielding, the core bounced and recollapsed in ~44 Myrs and ~111 
Myrs after the star has turned off in Phase I and Phase II, respectively (Table Hj 
Qualitatively, our calculation without H2 self-shielding agrees with the result of 



Abel et al. 



( 120061 ). that collapse in the target halo occurs after the source dies. We find that SIMF is 



the main mechanism for the formation of II2. Initially, the strong LW band photons destroy 
molecules in the core. As the shock propagates inward, however, boosted density and tem- 
perature of the post-shock gas enhances the molecule fraction (equation [26ll. and increases 
the H2 column density. As the shock front accelerates, SIMF occurs, and newly created H2 
is protected from the LW band photons because of increased self-shielding. If self-shielding 
is not accounted for, however, this H2 is destroyed and never restored, so collapse does not 
proceed during the lifetime of the source. 

We conclude, therefore, that neglecting H2 self-shielding in calculation explains why 



Abel et al 



(I2OO6I ) observes a delayed collapse. The quantitative disagreement between our 
collapse times (when we neglect self-shielding) and theirs may originate from the difference 
in the structure and chemical abundances of the target halo wh e n the source irradiates it. 



How do our results compare with those of 



Susa fc Umemural (120061 )? A fundamental dif- 



ference exists other than the fact that their work is limited to subclumps of a halo that 
hosts a Pop III star. They interpret the shock only as a carrier of negative feedback effect, 
while the shock, in our case, delivers both the positive and negative feedback effects. In 
their shock-driven evaporation (Model C) case, the collapsing core eventually fails to col- 
lapse, because the shock heats the core before it finishes collapse. Their successful collapse 
case (Model B) is simply an unaltered collapse: an already collapsing core finishes collapse 
before the shock front reaches the centre. On the other hand, we have observed expedited 
collapses as well as delayed or failed collapse. Such expedited collapses we observe are truly 
positive feedback effects. Quantitatively, because of their limited interpretation of the role 
of the shock, they argue that only regions with hydrogen number density nn ^ 10^~^cm~^, 
high enough to finish collapse before the shock front reaches the centre, can collapse under 
the influence of Pop III starlight. On the contrary, we find, for instance, that regions with 
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no radiation 


self-shielding 


no self-shielding 


Phase I 


11.2 


1.1 


47 


Phase II 


2.7 


1.3 


114 



Table 4. Collapse time (in units of Myrs) of a subclump with M = 5.5 X 10^ Mq irradiated by a Pop III star at distance 
D = 50pc (Fq = 600). For both Phase I and Phase II, we show how a case with a proper treatment of H2 self-shielding (2nd 
column) differs in collapse time from a case without self-shielding (3rd column) and a case without radiation. When H2 self- 
shielding is properly treated, collapse occurs in ~ 1 Myr, before the neighbouring Pop III turns off, while when H2 self-shielding 
is neg lected, collapse occurs after the star turns off, which is qualitatively consistent with the simulation results by .Abel et al.l 
l l2006r) . 



nn ~ 30 cm~^ - core density of TIS haloes in Phase I - can cool and collapse even after the 
shock front has reached the centre. As the shock-front accelerates and delivers strong positive 
feedback effects in the small core region, high resolution is required to pr oduce this mecha- 



Susa fc Umemura 



nism in simulations. The relatively poor resolution of SPH simulations by 
( 120061 ) might have prevented them from fully resolving the shock structure in the core, and 
potentially producing the positive feedback effects. 

Our result indicates that secondary star formation may occur even in subclumps of 
the host halo, which are subject to much stronger radiative feedback than isolated, nearby 
minihaloes. We have shown in this section that H2 self-shielding is important even at this 
high level of ionizing {Fq = 600) and dissociating (Flw = 2 x 10~^^ ergs~^ cm~^ Hz~^) fluxes. 
It is even more surprising because the collapse is expedited and coeval formation of Pop III 
stars in the same neighbourhood is possible. The naive expectation of negative feedback 
effect of a Pop III star in its neighbourhood, therefore, should be revisited. 



7 SUMMARY/DISCUSSION 

We have studied the radiative feedback effects of the first stars (i.e. Pop III stars) on their 
nearby minihaloes, by solving radiative transfer and hydrodynamics self-consistently using 
the 1-D spherical, radiation-hydrodynamics code we have developed. The results can be 
summarized as follows: 

• We identified the minimum collapse mass, namely the mass of minihaloes which are 
able to have a core which cools and collapses in the absence of external radiation. We find 
that Mc^min ~ 7 X 10'^ M0 at z = 20. In determining Mc^min, we applied two criteria. First, 
the collapsing region should reach rin = 10^ cm~^ to be considered as a collapse. Second, 
this should occur within the Hubble time. The minimum collapse mass we find roughly 
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agrees with that of iMachacek et al.l (120011 ). where the AMR scheme they used seems to have 



resolved the inner structure of minihaloes. 

• Minihaloes could have been in very different stages of their evolution when they were 
irradiated by a Pop 111 star. We used two different initial conditions to represent such phase 
differences. In Phase 1, chemical abundances have not yet evolved away from their IGM 
equilibrium values. This stage is characterized by low H2 fraction, ?/h2 ~ 2 x 10^^ and high 
electron fraction, x ~ 10~^ at the centre. Haloes can be irradiated in Phase II, which is the 
state of these haloes evolved from Phase 1, where x has dropped to 10~^ by recombination. 
Phase II is characterized by high H2 fraction ?/h2 ~ 10~^ — 10^'^, low electron fraction 
X = 10~^, and core density higher than that of Phase I. 

• Within our parameter space, the I-front is trapped before reaching the core in all cases. 
Ionized gas evaporates, and a shock-front develops ahead of the I-front and travels into 
the core. The shock front leads to both positive and negative feedback effects. A boost in 
density and temperature by a shock increases the H2 formation rate. In some cases, the shock 
accelerates and obtains a temperature above 10^ K, which is high enough to drive coUisional 
ionization, which then leads to a further boost in H2 fraction. The high temperature and 
kinetic energy delivered by the shock, on the other hand, tries to disrupt the gas. The nett 
effect is either 1) an expedited collapse, 2) delayed collapse, 3) neutral (unaffected) collapse, 
or 4) a disruption, depending upon the flux, halo mass, and the initial condition when 
irradiated. 

• At the moment of collapse, halo profiles under radiation are almost identical to those 
without radiation. Density profiles of different mass haloes are well fit by different power- 
law profiles, p oc r'"", where w = 2.5, 2.4, 2.3, and 2.2 for M = 10^ 2 x 10^ 4 x 10^ and 
8 X 10^ Mq, respectively. Some variation in temperature profile exists at the low- mass end, 
M = 10^ Mq, and the high-fiux end Fq = 46.3 {D = 180 pc). 

• Overall, the radiative feedback effect of Pop III stars is not as destructive as naively 
expected. Minihaloes with M > [1 — 2] x 10^ Mq are still able to form cooling and collapsing 
clouds at their centres even in the presence of radiation. A simple explanation is possible 
for such behaviour. In Phase I (low yn^ and high x), radiation can easily dissociate H2 
while the source is on, but after the source dies, high electron fraction allows H2 formation. 
On the contrary, in Phase II (high and low x), B.2 is more easily protected against 
the dissociating radiation because the higher H2 column density provides self-shielding and 
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compression increases the formation rate. The situation becomes more comphcated, however, 
by other feedback effects which will be described in the following bullets. 

• Within our parameter space, haloes that are irradiated at Phase I experience expedited 
collapse predominantly for 10^ < M/M© ^ 8 x 10^, except for the delayed or neutral collapses 
occurring at the low mass/high flux and the high mass/low flux extremes (e.g. for M = 
10^ Mq at Fq = 46.3 and for M = 8 x 10^ M© at Fq = [1.5, 5.14]). 

• Haloes that are irradiated at Phase II show a more complicated behaviour. In this 
case, unaffected collapse is more frequent, in general, at high and intermediate masses, 
while for M = 10^ M©, core collapse is now reversed at any Fq. Delayed collapse occurs for 
M = 2 X 10^ M© at Fo = 46.3. Unaffected collapse occurs for M = 8 x 10^ M© for any Fq, 
and for M = 4 x 10^ M© at Fo<11.6. Otherwise, for intermediate mass, collapse is either 
neutral or expedited. 

• We first find in this paper that coeval formation of Pop III stars is possible even under 
the influence of ionizing and dissociating radiation from a first star. This occurs either as 
an expedited collapse or an unaffected collapse. Among those parameters explored in this 
paper, expedited collapse occurs during the lifetime of the source star when a halo of mass 
M = 2 X 10^ M© in Phase I is irradiated by a Pop III star at a distance D = 180 pc 
{Fq = 46.3). Unaffected collapse occurs for haloes of mass M = 8 x 10^ in Phase II during 
the lifetime of the source star for all different distances (fluxes). 



• Extending our parameter space to include a specific case studied by lAbel et al.l (120061 ). 
a minihalo merging with a halo hosting a Pop III star, we find that the coeval formation 
of Pop III stars is possible even in this high i onizing (-Pn ~ 600) and dissociating (-Flw ~ 



2 X 10~^® erg s~^ cm~^ Hz~^) flux case. While lAbel et al.l (120061 ) find that the secondary 
star formation in this target halo occurs after the first star dies because of H2 destruction 
by photodissociation, we find that the minihalo core collapse is expedited to form a star 
in ~ 1 Myr, long before the first star dies, due to the SIMF and H2 self-shielding. This 
discrepancy comes from the fact that we account for the effect of H2 self-shielding, while 
they do not. A proper treatment of II2 self-shielding is important even for such a high flux 
regime, because the central H2 fraction can reach y^^ ^ 10~3 due to the SIMF and strong 
H2 self-shielding is possible due to newly created H2. 



We find the minimum collapse mass Mc^uim ~ 7 x 10^ M© at 2; = 20 wi thout radiation. 
While our result agrees roughly with that of the 3D AMR simulation by iMachacek et al. 
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2001) 



2000 



, discrepancy 



jecom es larger with those of 3D SPH simulation results (e.g. 



Yoshida et al 



( iTegmark et al 



Fuller fc Couchman 



20031 ) and a semi-analytical calculation using a uniform-sphere model 



19971 ). This implies that the central region of haloes should be resolved well 



in order to quantify the minimum collapse mass exactly. 

What does the result of our paper imply for the "first" H II region created by Pop III 
stars? Because a significant fraction of nearby minihaloes can host second generation stars 
within the first H II region, it is possible that such a subsequent star formation may at least 
keep the first H II regions ionized. It may even be possible that individual H II regions grow 
and overlap, thus finishing the first cosmological reionization. A semi-analytic calculation 
of minihalo clustering around high density peaks, for example, might allow us to quantify 
how fast and how big such bubbles can grow. Without secondary star formation, this would 



simply be a relic H II region in which gas recombines and coo^ 
possibly with metal enrichment from supernova explosion (e.g. 



s after the source s tar dies, 



Bromm et al. 



20031 ). 



We found that the minimum collapse mass is ~ 1 — 2 x 10^ Mq even in the p resence of Pop 



III star light. Such a low value may affect the reionization history significantly. 
2006a) estimates that the instantaneous ionized mass fraction at 2; = 20 is 



Alvarez et al. 



0.1, if indi- 
vidual ~ 10^ Mq haloes host one ~ 100 Mq Pop III star each. If the typical mass scale of 
host haloes is ~ 10^ Mq instead, as the number density of haloes would be roughly 10 times 
as big as that for M ~ 10^ Mq, Pop III stars alone would be able to finish cosmological 



reionization at z 



2C 



_ New reionization sources will form later in more massive haloes 
with Tvir ^ 10^ K, which will host a region cooling by the hydrogen atomic cooling. Depend- 
ing upon how fast such transition occurs, the global reionization history will have different 
characteristics (e.g. monotonic growth of ionization fraction vs. double reionization). 

In this paper, we have considered only the radiative feedback effect. Pop III stars, how- 
ever, may exert additional feedback effects. The H II region developed by a Pop III star 
inside the host halo breaks out as a "champagne flow" inside the host minihalo, where the 
I-front separates from the shock-front and runs ahead, transforming from D-type to R-type. 
The shock front left behind also expands into the IGM and nearby minihaloes would be 
encountered by this shock-front ultimately. Other feedback effects will come from supernova 



^ This argument is based upon the fact that the comoving number density of haloes, M dn/dM, is roughly proportional to 
. The minihalo population, however, might hav e been severely re duced by the "Jeans-mass filtering" inside ionized bubbles 



created around rare, but more massive objects (e.g. 
negligible contribution to cosmic reionization. 



Iliev et al 



20061 ) . in which case sources hosted by minihaloes would make 
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explosions. If the first star dies and explodes as a supernova, both dynamical and chemical 
feedback effects would alter the fate of nearby minihaloes, as well. 



How would the additional presence of H2 dissociating background radiation affect our 
results? In this paper, we have considered the effect of the radiation from an individual 
nearby Pop III star, whose SED takes a black body form for a short lifetime (~ 2.5 Myrs). 
This is the case appropriate to the earliest star formation. It is valid whenever a minihalo 
resides in a place and time where the background from other, more distant stars is negligible. 
On average, however, the mean free path to H2 dissociating radiation is greater than that 
for ionizing radiation prior to reionization, so the situation can arise in which the ionizing 
radiation from distant sources is filtered out but the UV radiation in the LW bands is not. 
Suppose a minihalo is under the influence of both Pop III starlight from a nearby star 
and a persistent background radiation field in the LW bands. In the absence of the nearby 
star, the dissociating background can only hinder the formation of II2 and its cooling. As 
such, the H2 fraction inside the minihalo when the nearby Pop III star starts to irradiate 
it would be lower than it would have been without the background. In this case, even if 
the background were intense enough on its own to prevent the minihalo from cooling and 
collapsing, the minihalo could still host a cooling core if H2 formed by the positive feedback 
from the Pop III star, despite the presence of the background. Indeed, this could occur 
frequently, because we find that a high electron fraction - and, thus a high II2 fraction - 
can be achieved by collisional ionization in the postshock region in many cases (SIMF; see 
Section [6.3. ip . This newly created II2 will then be easily protected from the dissociating 
background by self-shielding, since our simulation results show that this SIMF H2 survives 
even the much larger - albeit short-lived - flux of H2 dissociating radiation from a nearby 
star in our most extreme case, Flw ~ 2000 x lO^^^ergs"^ cm~^ Hz~^, as has been shown in 
Section 16.61 Thus, the background would then only prevent those haloes that cannot "host" 
this SIMF mechanism from cooling and forming stars. We will address this issue further in 
the future. 

As the focus of our paper is the fate of neutral cores of target haloes, in which the 
ionized fraction never exceeds ~ 10~^, we neglected processes which are relevant only when 
gas achieves high ionized fraction, such as H P cooling and charge exchange between He+ 



(He) and H (H+) (see e.g. lYoshida et al.ll2006l ). These processes may be important, however. 



in the relic H II region outside the target minihalos. For instance, HD cooling may cool gas 
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down below the H2 cooling temperature pla teau, Th^ ~ 100 K, i 



start from a highly ionized initial state (e.g. 



Johnson fc Bromm 



H9. f ormation and cooling 



20061). 



We chose two different evolutionary phases of nearby minihaloes as our initial conditions. 
A more natural way to address this problem is to use the structure and chemical composition 
of minihaloes and IGM from 3-D, chemistry- hydrodynamics calculation. We intend to ex- 
tend our study in a more consistent manner by combining a 3-D, chemistry-hydrodynamics 
simulation and the 1-D, radiation- hydrodynamics simulation in the future. In this paper, 
we simply adopted a model for virialized haloes (TIS pr ofile). In the future, we will also 



implement a more realistic growth history of haloes (e.g. 
for the dynamical effect of mass accretion. 



Wechsler et al. 



20021 ) to account 
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APPENDIX A: NUMERICAL METHOD AND CODE TESTS 

Here we describe the finite-difference scheme used for our 1-D spherical, radiation-hydrodynamics 
code. The subscript, unless noted otherwise, denotes the position of a shell. The superscript 
denotes the time. For instance, Pj+1/2 is the zone-centred density of shell j + 1 at time 
and rj is the zone-edge-centred radius of shell j at time t". 



Al The Gas Dynamical Conservation Equations 

Hydrodynamic conservation equations for the baryonic component ( eqs. [3] - [5]) are solved 
following the finite-difference scheme by iThoul fc Weinberg] (119951 ). We first update the 
velocity and position using the so-called "leap-frog" scheme, so that the velocity and the 
position are staggered in time: 



n+1/2 n-1/2 



and 



Pj-1/2 



dm a 



+ 



(Al) 
(A2) 
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which are second-order accurate. As the mass of each shell is conserved for such a Lagrangian 
scheme, density is updated following 

= dmj+i/2 

^^■+1/2 (4/3)7r[(r;+i)3 - (r^+y] ' ^""^^ 
In these equations, 

= ^ ^^"+1/2)^ (A4) 

and 

drrij = ■^(rfmj_i/2 + dmj+1/2). (A5) 
We then advance the energy by 

n+1 _ n , _ „n / _J: L 

^j+1/2 — W+1/2 ^+1/2 n+1 

\ri+l/2 Pj+1/2^ 

+i£Z^Wrfr+i/2. (A6) 
Pi+1/2 

Shocks are treated with the usual artificial viscosity technique. The pressure in the 
momentum and energy conservation equations is replaced hj P = p + q, where 



,n+l ^ 



n+1/2 n+1/2 



'^i/pr+\%-i/pr+i/2 

^iv:T-vr'% (A7) 

if vl^i^"^ — f "^^^^ < 0, and g = otherwise. We use Cq = 4, which spreads the shock fronts 
over four or five cells. 

Dark matter shells are also updated according to equations flAip - flA7p - note that we 
use fluid approximation as described in Section [2T2] -. except that the heating/cooling term is 
zero in equation (]A6|) . Note that the dark matter shells are allowed to have effective shock in 
our fluid approximation, and therefore we need to compute the artificial viscosity when dark 
matter shells are converging (equation IA7I1 . as in the case of the baryonic gas component. 

A2 Time Steps 

Time step for the finite-differencing is chosen such that important fluid variables do not 
change abruptly. The relevant time scales are the dynamical, sound-crossing (Courant), 
cooling(heating), and species-change time scales. In addition, to ensure that the fluid shells 
do not cross, we also adopt a shell-crossing time. 

dt = mm{dtdyn, dtcom, dtcool, dtspec, dtyei} 
© 0000 RAS, MNRAS 000, 000-000 




mm < 



(A9) 

(AlO) 

(All) 
(A12) 
(A13) 

where q, cq, Cc, Cgp, and c„ are coefficients that ensure accurate calculation of the finite 
difference equations. We use Cd = 0.1, cc = 0.1, Cc = 0.1, Csp = 0.1, and Cy = 0.05. 

In practice, we frequently find that dt^y^ can be very small compared to other time scales. 
We sometimes disregard dtdyn in order to achieve computational efficiency. We confirmed, 
especially in our problem, that such a treatment does not produce any significant discrepancy 
from a calculation with dt^yn considered. When the virial temperature of a halo is close to 
the cooling temperature plateau, for instance, dt^yn must be irrelevant because gas would 
be almost hydrostatic. 



A3 Radiative Transfer 

For the radiation field generated from a point source at the centre, the radiative rate coeffi- 
cient of species i at radius r is given by equation ([T5l) . Finite-differencing this rate coefficient, 
however, requires some caution. For the baryonic shell at position j (smaller j means closer 
to the centre) whose inner edge and outer edge have radii ?"j-i/2 and rj+1/2, respectively, the 
incident differential fiux at the outer edge is F"*(rj+i/2), and one could naively calculate 
the rate coefficient of species i by 

Hr,)=j^ du . (A14) 

As mentioned already in Section 12.3.21 and Section 12.4.11 however, this expression may 
not yield an accurate result when the shell k is optically thick. In this case, F^, may change 
substantially over the shell width, and equation (lA14p might overpredict the ionization 
rate by applying a constant fiux over the shell width {Arj = Tj+1/2 — Tj_i/2). One may, in 
principle, choose to set up the initial condition such that all shells are optically thin. However, 
such a scheme can be very expensive computationally, especially when collapsed haloes 
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are treated. In ord e r to r esolv e this prob l em, w e use the "photon-conserving scheme" by 



Razoumov fc Scott 



(119991 ) and 



Abel et al. 



(119991 ). In this treatment, the number of photons 
that are absorbed in a shell is the same as the number of ionization events. Equation (1A14I1 
can then be re-written as 



ki{rj) 



du- 



du- 



hu 



i shell J 



hv 



1 _ Q~^n,^{rj) 
HiAVj 



(A15) 



where L^^*(r) = 47rr^Fj^^*(r), Arj,^(rj) = UiArjai^i, is the optical depth of a shell k on 
a species i, and Kheiij — 47rr|Arj is the volume of the shell. Note that when Ar^ ^ 1, 
equation (IA15P becomes equivalent to equation (IA14p . For each species, the corresponding 
radiative reaction rate is calculated by quadrature, by summing the integrand in equation 
(]A15p . then summing over the frequency to obtain the nett radiative reaction rate. 



A4 Nonequilibrium Chemistry 



As described in Section 12.51 in order t o upd ate the abundance of species i, we adopt the 



finite difference scheme by 
updated by 



Abel et al. 



n+l 



( 1l997l ). Based upon equation (fT7|) . each species i is 



(A16) 



Abel et al. 



l + D^+\T,{nj})dt^+y^ ' 
where the species {%} is the previously updated value in the order given by 
( 19971 ) (note that the letter n (n + 1/2, n + l) in superscript denotes the time (t""'"^/^, 
fn+iy rjj^g order they find to be optimal is H, 11+ , He, He"*", He"'""'" and e~, followed by the 
algebraic equilibrium expressions for H~ and H"'", and finally H2, again by equation ( ]A16p . 



A5 Numerical resolution 

In practice, we use 500 dark matter and 1000 fluid shells sampled uniformly (in radius) 
from the centre to the truncation radius rtr. We put a small reflecting core at the centre 
with negligible size, namely rcore = lO'^^rtr. Such a core is found to be useful in reducing 
undesirable numerical instability at the centre. Our choice is conservative enough not to 
affect the overall answer. 

A wide range of radiation frequency (energy), /iz/ ~ [0.7 — 7000] eV, is covered by 100, 
logarithmically spaced bins, AE/E ^ 0.04, together with additional, linearly-spaced bins 
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Izotov & Kolesnik f 1984) 


16 


H- 


f H+ 


2H 


Dulev & Williams (1984) 


17 


H- 


f H+ 


^H++e- 
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Schneider et al. (1994) 
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Osterbrock (1989) 
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de Jone (1972): Shaoiro & Kane (1987) 
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Dunn (1968) 
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Oneil & Reinhardt ( 1978) 
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f 7- 


2H+ + e- 


Bates & Opik (1968) 


27 


Ha + 7 ^ 
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Sectionl2.3.1l Draine & Bertoldi (1996) 



Table Bl. Reactions and the corresponding references. 

where radiative cross sections change rapidly as frequency changes. About a dozen hnearly 
spaced bins at each of those rapidly changing points turned out to produce reliable results. 



APPENDIX B: RATE COEFFICIENTS 



In Table IBll we list the chemical reaction rates we implemented in our code and the cor- 
responding references. The rate coeffici e nts ( 1-19) and radiative cross sections (20-26) are 
mostly from the fit by lShapiro fc Kangl (Il987l ). except for a few updates. 



APPENDIX C: CODE TESTS 



We now extend the description of our code test problems in Sect ionl2.6land show t 



le results. 



(A) The self-similar, spherical, cosmological infall problem (iBertschingerl Il985l ) : A point 



mass, if placed in an unperturbed Einstein-de Sitter universe, will make all particles around 
it to be gravitationally bound, leading to a successive turnaround and collapse of spherically 
shells. Infalling matter will be shocked and form a virialized structure, whose profiles are 
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well described by a self-similar solution. We restrict ourselves to purely baryonic fluid with 
the ratio of specific heats 7 = 5/3. 

The turnaround radius rta, at which the Lagrangian proper velocity of a shell is zero, 
evolves as 

rut)={^y''\s^R^y''mf', (ci) 

where 5iRf defines the seed mass Sm added to the Einstein-de Sitter universe, 
4 

Sm = -7ipu,ARf, (C2) 

where the initial cosmic mean density pH,j = l/(67rGtf) at t = t,. The shock radius is 
a constant fraction of rta: r^it) = 0.338976 rta(t) for 7 = 5/3. The dimensionless radius 
A = r/rta and the dimensionless density D = p/pu, where the cosmic mean density pn = 
l/iQirGt^), satisfy the unique Bertschinger solution. In Fig. ICll we show the density profiles 
and rs{t), obtained from the simulation with 6iRf = 1.84 x 10^^ cm^, tj = 5.572 x 10^^ s. 

( B) The self- similar blast wave from a strong, adiabatic point explosion in a uniform 
gas (jSedovl Il959l ): A point explosion drives a self-similar blast wave through the initially 
static, uniform medium. A strong shock is generated, and rs(t) = ^0 (— ) ^ t"^^^, where E is 



the thermal energy of explosion, po is the initial density, and ^0 is a dimensionless constant 
determined by 7. For 7 = 5/3, ^0 = 1-152. We use E = 1.053 x 10^"*^ erg, 7 = 5/3, and 
Po = 2.5626 X lO'^^'cm"^ for simulation results displayed in Fig. ICll 

(C) The propagation of an I-front from a steady point-source in a uniform, static medium: 
This is the case where the classical description of the Stromgren radius is plausible, since 
gas is forced to remain static, and photoionization and recombination are the only physical 
processes determining the ionized fraction. The I-front from a point source with A^,, number 
of ionizing photons evolves as 

ri(t) = Rs{l- exp(-t/trec))'/' , (C3) 

1 /3 

where Rs = [3A^*/(47rn^a)] is the Stromgren radius, tree = I/I'^h^) is the recombination 
time, and a is the recombination rate coefficient. We adopt A^* = 10^^ s~^, = 10cm~^, 
and a = 1.05 x 10~^^ cm'^ s~^. For this test, we use a monochromatic light whose frequency 
is slightly above the hydrogen ionization threshold. 

( p) the gas-d ynamical expansion of an H II region from a point source in a uniform 
gas ( Lasker 19661 ): The I-front, initially propagating as a weak R-type front into a uniform 
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Figure CI. Code test res ults. From top to bottom, simulation results of (A) the self-similar, spherical, cosmologica l infall 
probl em ijBertsch ingcr' 1985), (B) the self-similar blast wave from a strong, adiabatic point explosion in a uniform gas (Sedovl 
Il959h . (C) the propagation of an I-front from a steady point-s ource in a u niform, static medium, (D) the gas-dynamical 
expansion of an H II region from a point source in a uniform gas l|Laskedll966h . and (E) the gas-dynamical expansion-phase of 
the H II region from a point-source in a nonun iform gas whose density varies with distance r from the source as r~™, w = 3/2 
ijPranco. Tenorio-Tagle. &: Bodenheimei]ll99Ci ) are displayed, respectively, as described in text. In each row, some early (circle) 
and late (square) snapshots of density profiles (neutral fraction profile in case C) are shown in the left panel, while the evolution 
(cross) of shock radius (r^) or I-front (rj) are shown in the right panel. Data points (circle, square, cross) are compared with 
the analytical prediction (solid line), in case analytical solutions exist. Note that for test (D), ri is plotted against i = t — tc, 
where tc is the time when dri/dt = cj. In this case, the analytical solution for rj is valid only for i > 0. On left panels, snapshots 
are shown at t = 5.5 X lO^'^s and 1.1 X for (A), t = 5 X lO^'^s and lO^^s for (B), t = 1.4 X 10''yr and 6.1 X lO^yr for (C), 

t = 9x lO^yi and 2 X lO^^yr for (D), and t = 4.3 X Iff^yr and 8.6 X 10**yr for (E). 



medium, slows down and travels as a D-type front, developing a shock front ahead of it. The 
I-front evolves as 

4/7 



ri(t) = Rs,i 1 




(C4) 

where Rs,i is the initial Stromgren radius, t = t — tc is the time measured from the moment 
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tr. when dr]/d t = ci, and ci = (p/pY^'^ is the isothermal sound speed of t he ioni z ed ga s 
( Spitzer 1978 ). We adopt A^, = 2.45 x lO^^s-^ and = 6.4cm-3. Following llaster ( 1966 ). 
we force temperature of the ionized gas to be 10'^ K, which gives ci = 12.86 km/s. 

(E) The gas-dynamical expansion-phase of the H II region from a point-source in a 
nonuniform gas whose density varies with d istance r from the source as r""", w = 3/2 



199( 



1^: This case is similar to the case (D), except 
that the density follows a power law, tih oc r~^. Inside the core radius rc, the density is 



( iFranco. Tenorio-Tagle. fc Bodenheimer 



constant at nH,c- The I-front evolves as 
7-2w / 12 \ dt 



ri{t) = K 



1 + 



4 



9-4W 



(C5) 



where Rw is the size of the initial H II region obtained by equating the ionization rate and 
the recombination rate. For instance, when w = 3/2, 



R 



3/2 = r^exp 



1 



(C6) 



1/3 



where Rs = / {Ann^^^a) 

If w < 3/2, the shock front always travels ahead of the I-front. If w > 3/2, however, 
the shock front is overtaken by the I-front, which soon runs to infinity in this "champagne" 
phase. We restrict ourselves to this critical exponent w = 3/2. From equation ( 1C5I) . we obtain 
ri{t) = R3/2 (1 + 2C1/-R3/2). In our simulation, we use nH,c = 2x 10^ cm" 



2.1 xlO^'' cm. 



5 X 10^^ s 1, and a = 2.6 x 10 ^'^cm'^s ^. Temperature of the ionized gas is set at 



'13 



T = 8000 K, such that cj = 11.5 km/s. 
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